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Abstract: We perform a global Bayesian fit of the phenomenological minimal su- 
persymmetric standard model (pMSSM) to current indirect collider and dark matter 
data. The pMSSM contains the most relevant 25 weak-scale MSSM parameters, 
which are simultaneously fit using 'nested sampling' Monte Carlo techniques in more 
than 15 years of CPU time. We calculate the Bayesian evidence for the pMSSM 
and constrain its parameters and observables in the context of two widely different, 
but reasonable, priors to determine which inferences are robust. We make infer- 
ences about sparticle masses, the sign of the n parameter, the amount of fine tuning, 
dark matter properties and the prospects for direct dark matter detection without 
assuming a restrictive high-scale supersymmetry breaking model. We find the in- 
ferred lightest CP-even Higgs boson mass as an example of an approximately prior 
independent observable. This analysis constitutes the first statistically convergent 
pMSSM global fit to all current data. 
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1. Introduction 



There is currently an expectation that with the start of the Large Hadron Colhder 
(LHC), high energy physics will soon enter a new phase highly dominated by new 
data that could imply physics beyond the Standard Model (SM). Over the past years, 
low energy supersymmetry (SUSY) has become the standard approach to study the 
potential physics beyond the SM, mostly because of its natural power to address the 
hierarchy problem. The 124 Lagrangian parameters of the minimal supersymmetric 
extension of the SM (for a recent review, see [|I|) makes its phenomenological study 
impractical. It may well be that the mechanism that mediates SUSY breaking to 
the observable sector provides relations between many of these parameters. Unfor- 
tunately, however, there are many different mediation mechanisms in the literature, 
with no one clearly preferred. 

In order to extract computable information, many works have reduced the num- 
ber of parameters by truncating to a handful of soft-breaking parameters at a high 
energy scale. The remaining set of parameters are used as boundary conditions for 
renormalisation group equations (RGE), which are run down to the weak-scale. A 
large amount of minimal supersymmetric standard model (MSSM) based studies 
have been carried out in the way described above. Most of them were performed 
in the context of the constrained MSSM (CMSSM, also sometimes called mSUGRA 
for minimal supergravity) set-up which have only 4 independent non-SM parameters 
(and a sign choice). Many groups have been pursuing a programme to fit this model 
and identify regions in parameter space that might be of interest with the forth- 
coming LHC data. See for instance [0, |^, ^ |], |^, |^. Complete scans over up to 8 
free parameters of CMSSM with a combined treatment of likelihoods from different 
experimental constraints were possible with Markov Chain Monte Carlo (MCMC) 
sampling techniques P, ||, |, |10|, |ll], |12|, However, the truncation to a handful 



of parameters in the CMSSM, is at best a very strong assumption, and most likely 
over-restrictive. 

There are two directions that can be followed to properly study low-energy super- 
symmetric models. The top-down approach has been tried over the years by deriving 
the otherwise free parameters from an ultraviolet extension of the MSSM. Models of 
unification, different sources of SUSY breaking, such as gravity and gauge mediation 



and classes of string compactifications [|^, |T6], [T^, |T^, 0, |0], ^ |2^, |j have been 
used to provide high energy expressions for the soft breaking parameters. Soft SUSY 
breaking terms can be computed at energies as high as the GUT scale of ~ 10^^ GeV, 
and renormalisation group running to the TeV scale allows contact to be made with 
potential quantities of interest such as sparticle masses. Recent progress in moduli 
stabilisation in string theory has made this approach more concrete and calculable 
with explicit results for some classes of models. This is very encouraging but usu- 
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ally the string derived models fill only a small subset of the full CMSSM parameter 
space, which could make them impossible to differentiate from the CMSSM. The 
proliferation of SUSY breaking mechanism set-ups mean that analyses where only 
one is picked tend to be very specific, with a rather limited range of applicability. It 
is desirable to side-step such extreme model dependence with a different approach. 

Alternatively, one can use a bottom-up approach to low-energy SUSY. In this case 
the soft breaking parameters are considered at the SUSY scale without referring to 
their high energy origin. All the parameters can in principle be varied over the ex- 
perimentally allowed range and compared with potential observations at the LHC 
and other experiments. It is a formidable task to consider all the 124 parameters, 
due mostly to computing limitations. An interesting compromise is the pMSSM p4| 



which we consider here. In this model, the number of free parameters is 20 soft- 
breaking parameters (and a ±1 parameter) plus 5 SM ones and are selected by the 
requirements of consistency with unobserved flavour changing and CP violating pro- 
cesses. However, even this simplified version of the MSSM requires a lot of computer 
power to be analysed in and details. The ability of future collider measurements to 



constrain the pMSSM has been estimated using MCMC methods in Refs. |]25|, ^ 



At present, SUSY forecasts for the LHC necessarily contain large uncertainties. In 
particular, there is a strong model dependence on the mechanism for SUSY break- 
ing. Realistic predictions need guidance from direct and precise (collider and other 
related) experimental data. Interestingly, the converse is also very important: the 
experiments need unbiased phenomenological guidance about the expected nature or 
properties of SUSY. This is what we aim to accomplish, eventually. 

The purpose of this article is to perform a global fit of the pMSSM and make SUSY 
forecasts for collider and dark matter search experiments using Bayesian statistics 
methods. For given prior probability and likelihood densities, Bayes theorem provides 
the way to extract the posterior probability density for the parameters. It can also 
be used for model comparison when enough data is available. This formalism has 
been used in many fields of science, including cosmology (see |]27| , p8| for recent 
reviews). Markov Chain Monte Carlo (MCMC) and related techniques have recently 
begun to be used to perform Bayesian inference on supersymmetric models. The 
increasing access to large scale computing power and improved methods of calculation 
are making these techniques more manageable with time and we have been able to 
tackle the relevant parameters of the pMSSM. 

The complete and simultaneous scan of the 25 parameters and a sign^ for the 



pMSSM construction was performed using the MultiNest program p9|, pO|. At 



^Where we refer to 25 parameters, we shall really mean 25 continuously varying parameters plus 
the one discrete sign choice. 
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the heart of the algorithm is the nested samphng technique |31] that revolutionised 
computational Bayesian inference by prioritising a computation of the Bayesian evi- 
dence rather than solely on computing the posterior probability distribution function 
(PDF) of model parameters (although the latter is obtained at no additional cost), 
as is the case in traditional Monte Carlo algorithms (e.g. MCMC). We will review 
the nested sampling method in an appendix. In simple terms, similar to the MCMC, 
nested sampling is an iterative Monte Carlo method that, starting with a relatively 
small number of points (a few thousand in our case), it produces a list of a large 
number of points (10^ — 10*^ in our case) ordered in increasing likelihood. We used 
it here because it computes both the evidence and the posterior PDF, as opposed 
to traditional MCMC methods which only calculate the posterior PDFs. Further 
the MultiNest algorithm is efficient in handling complicated problems with multi- 
modal/degenerate posterior distributions. 



We emphasise that the current situation with no direct sparticle measurement 
data yet from LHC makes the issue of prior dependence critical. For this reason, 
it is expected that extraction of prior independent information of our analysis will 
be difficult. Interestingly, however, we find some results with approximate prior 
independence. In order to illustrate the issue of prior dependence of results we 
consider priors that are flat in the parameters themselves ('linear') and flat in the 
logarithm of the parameters ('log') priors. As usual in Bayesian statistics prior 
dependence should not be understood as a drawback but as a positive feature that 
can be used to determine when enough data is available to unambiguously make 
inferences. It is expected that, if SUSY is discovered, the addition of LHC sparticle 
mass data will relax any prior dependence, and so an analysis along the lines of ours 
could be used to extract prior independent information. 



The amount of work related to this project required very efficient algorithms and 
access to high performance computing. We used the University of Cambridge super- 
computers: COSMOS from the Department of Applied Mathematics and Theoretical 
Physics (DAMTP) and the Darwin cluster from the High Performance Computing 
Service (HPC). The final run was made in terms of 60 12-hour jobs, each correspond- 
ing to a cluster of 128 CPUs on HPC and 40 8-hours jobs, each corresponding to 
a cluster of 64 CPUs, on COSMOS (making a total of more than 15-year standard 
CPU time). Some results of a complete and independent (from the one we present 
here) runs with fewer experimental constraints were presented by one of us (SSA) 
at the SUSY 2008 conference and reported in Ref. ||3^. While we were upgrading 
our analysis, a study of randomly scanned pMSSM points appeared similar in 
philosophy to Ref. 0]. In the conclusions, we contrast the aims and methodologies 
of our work with these. 
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In order to make a self-contained presentation we briefly describe Bayesian infer- 
ence and relevant terminologies in Section We construct the elements needed for 
inference in the context of the pMSSM in Section ^. The experimental constraints or 
observables used are described in Subsection |3^ . The sampling procedure and the 
different high energy physics software used to predict the observables are presented in 
Subsection In Sections ^ and ^ we analyse our results and then conclude. In the 
Appendix we briefly review the nested sampling method and how the MultiNest 
program works. 



2. Bayesian Inference 

Bayesian inference fits and plays an important role in the scientific process of data 
collection and modelling. It particularly deals with the steps that involve model 
fitting to data and the technique of assigning preferences to alternative models (model 
comparison). This subject is very important especially with the imminent start of 
the LHC experiments. Here we will give a short review of the basics of Bayesian 
statistics that are useful in our work. 

2.1 Bayes' theorem 

Consider a given model or hypothesis H (we shall take H to be the pMSSM) 
defined by some set of parameters (in our case, 25 parameters) 0. We wish to know 
the PDF P(0|D, H) of the parameters given the data D and the model set-up, H . 
P(0|D,if), being the parameters' PDF after confrontation with data, is called the 
posterior PDF. The likelihood, P(D|0, H) = L{@) is a measure of how well a model 
point predicts data set D. In order to calculate the posterior from the likelihood, 
one must assign some prior PDF to the parameters P{@\H) = 7r(0) to parametrise 
our uncertainty in them before the model is confronted with data. Bayes' theorem 
then describes how one may obtain the posterior from the other two PDFs and a 
normalisation constant P(D\H) = Z, the Bayesian evidence for the model in light 
of the data: 

The Bayesian evidence is given by 

Z = j L(0)7r(0)c/^0. (2.2) 

Here is the dimensionality of the parameter space; N = 25 for the pMSSM. Since 
the Bayesian evidence does not depend on the parameter values 0, it is usually 
ignored in parameter estimation problems and posterior inferences are obtained by 
exploring the un-normalised posterior using standard MCMC sampling methods. 
However, the evidence plays a central role in our discussion. 
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A useful feature of Bayesian parameter estimation is that one can easily obtain 
the posterior PDF of any function, /, of the model parameters 0. Since, 

P(/|D) = Jpif,&\B)d& = y"p(/|0,D)P(0|D)rf0 

= J (5(/(0)-/)P(0|D)d0 (2.3) 

where the probability chain rule is employed for the second equality and S is the 
Dirac delta function. Thus one simply needs to compute /(0) for every Monte 
Carlo sample and the resulting sample will be drawn from P(/|D). We make use of 
this feature in Sections ^ and ^ where we present the posterior probability PDFs of 
various observables used in the analysis of the pMSSM. 

Before proceeding to discuss prior distributions, we first briefly address ^, the dif- 
ference between the Bayesian and the frequentist approaches to inference. Bayesian 
inference (using Bayes' theorem) is a robust technique for updating prior knowledge 
or belief based on new data. It is unlike the frequentist approach where observations 
are viewed as random draws from some pool of possible observations such that the 
probability of an observation is the frequency observed with a large number of re- 
peated measurements. Frequentists usually focus on the likelihood and argue that 
the Bayesian approach is too subjective because of the use of priors. Bayesians reply 
that when many examine frequentist statistics, they are actually implicitly using pri- 
ors anyway. If one thinks of a region of parameter space with the lowest chi-squared 
values as being more likely than a region of parameter space with much higher chi- 
squared values, one is implicitly using some vague prior. The two approaches are 
asking different questions. Bayesians ask, "how likely is a given parameter value 
given the data?" , while frequentists ask, "how probable is the data, given certain pa- 
rameters?" In situations where the data is very informative the two approaches give 
the same results. We are interested in PDFs of parameters and so we use Bayesian 
statistics. By comparing results from different but reasonable priors, we obtain an 
estimate of how robust an inference is given current data. 

2.2 Priors 

The prior probability of the model parameters is a PDF that gives a subjective 
measure of our initial knowledge/ignorance about the values of the parameters of 
the model before the data are taken. Symmetries and physical observations or ex- 
pectations are usually a good guide to which priors to take. Information theory also 
provides a way of selecting priors by favouring those that maximise the entropy of 

^We follow the description in: 

http : //www. dsg.port . ac.uk/ valiviitaj/Lectures2006/CrittendeiiCMB2003 .pdf 
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the distributions. Some commonly used prior PDFs P{@\H) are: 



P{@\H) oc constant - the hnear prior, flat in (2-4) 
P{@\H) fx ^ - the Jeffreys prior, flat in log(0) (2.5) 
P{&\H) oc e-(®-®)'/2a2 _ ^j^g Gaussian prior. (2.6) 

The linear priors are often used for translational (like time and location) parameters, 
where there is no information to suggest that one value is preferred over others. The 
Jeffreys prior, also referred to as logarithmic (log) prior, is uniform in the logarithm 
of the parameter. These two priors are improper since they diverge when integrated 
over an infinite range. Our log and linear priors will be bounded by the requirements 
of perturbativity of the model, by passing previous direct sparticle search constraints 
and by the requirement of not too large fine-tuning in the Higgs potential parameters. 
These three criteria are sufficient to bound all 20 non-SM input parameters to a finite 
range. The Gaussian prior, on the other hand, is proper and integrable but requires 
previous experimental knowledge on a and 0. This is indeed the case for our 5 SM 
input parameters, and we use Gaussian priors for them (see Section I 



Assuming the parameters are independent, the resultant prior is obtained by the 
product of all the prior probability densities for each of the individual parameters. 
For instance, in the case of pMSSM with 25 parameters, 61,62, ... , 625, 



25 



P{@\H) = n{@) = l[n{6i). (2.7) 



1=1 



For our construction and analysis we are going to choose a linear prior measure 



for the pMSSM parameters described in Section This is because there are no 
observational evidence that hint to giving preference of some parameter region over 
others. We are going to determine bounds on the parameters from the fact that 
the parameter values have to be not far away from the TeV-scale in order avoid 
the little hierarchy problem. We are going to call such case/scenario as the linear 
prior. We then check the dependence of the results in our analysis on prior change 
by performing another analysis with a log prior. When we refer to a log prior, in fact 
the Jeffreys prior is only used for all parameters which only have positive bounds. 
Parameters which are allowed to take either sign present a problem with Jeffreys 
priors since the prior diverges at the origin. Therefore, the priors of such parameters 
are always taken to be linear. 

For parameter estimation, the priors become irrelevant once the data employed are 
powerful enough. This is already in evidence by comparing Bayesian CMSSM fits 



with those in similar supersymmetric models which have a lower number of free 
parameters, for example the large volume string compactification (LVS) scenario 
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The LVS scenario has two less free parameters than the CMSSM and current indirect 
data are aheady enough to make the result approximately prior independent. We 
may expect the addition of two precise, constraining, non-degenerate measurements 
(such as sparticle mass measurements from the LHC) to have the same effect upon 
the CMSSM. 

For model comparison, the dependence on priors always remains (although with 
more informative data the degree of dependence on the priors is expected to decrease. 



see e.g. [^j). Indeed this explicit dependence on priors is one of the most attractive 
features of Bayesian model selection. Refs. [|^, ^ identified prior distributions in 
high-scale CMSSM Lagrangian parameters. In particular, a Jacobian was defined 
to transform between derived parameters (such as tan/3) and more fundamental 
Lagrangian parameters from which they are derived. It is not our purpose here 
to find the 'most natural' prior because any such choice is necessarily subjective. 
Instead, we shall check the robustness of any inference under a reasonable variation 
of the priors. Such a check is especially required in model comparison hypothesis 
tests, which may be particularly sensitive to the particular choice of prior and its 
associated metric in parameter space 

2.3 Model comparison 

In order to evaluate and rank two alternative models Hq and Hi in the light of data 
D one needs to compare their respective posterior probabilities given the observed 
data set D, as follows^: use Bayes' theorem to relate the plausibility of Hi given 
the data, P(ifi|D), to the predictions made by the model about the data, P(D|ifi), 
and the prior plausibility of Hi, P{Hi). With this procedure one could construct the 
following probability ratio 

P{Hi\B) _ PiB\Hi)PiHi) _ ZiPjHi) 

P(/7o|D) P{Ti\H^)P{H^) ZoP{Ho)' ^ ' ' 

P{Hi) / P{Hq) is the prior probability ratio for the two models, it measures how much 
our initial beliefs favour Hi over H2. It is often set to unity but may occasionally 
require further consideration. The other ratio, P(D\Hi)/P(D\Hf)) measure how well 
the observed data were predicted by Hi and Hq. It can be seen from |2]^ that the 
Bayesian evidence takes the centre stage in Bayesian model comparison. As the 
average of likelihood over the prior, the Bayesian evidence is higher for a model 
if more of its parameter space is likely over some region of significant integrated 
prior ('prior mass') and smaller for a model with highly peaked likelihood but has 
a large prior mass associated with low likelihood values. Hence, Bayesian model 
selection automatically implements Occam's razor, since the prior PDF decreases 
for models with higher numbers of parameters. A more complicated theory will 



^Here we follow the description by 
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only have a higher evidence if it is significantly better at explaining the data than a 
simpler theory with less parameters. This technique was applied in to compare 
two CMSSM models: sign{fi) > versus sign{fi) < 0. We shall perform a similar 
comparison for the pMSSM in Subsection ^4.4| . The comparison of different GUT 
scale SUSY breaking models is also interesting . 



Another perspective on model comparison is in quantifying the consistency be- 
tween two or more data sets or constraints Different experimental con- 
straints may 'pull' the model parameters to different directions and consequently 
favour different regions of the parameter space. Any obvious conflicts between the 
observables are likely to be noticed by the 'chi by eye' method commonly employed 
to-date, but it is imperative for forthcoming high-quality data to have a method 
that can quantify such discrepancies. The simplest method for analysing different 
constraints on a particular model is to assume that all constraints and data provide 
information on the same set of parameter values. This can be considered as one 
hypothesis or model, Hi. This is the assumption which underlies the joint analysis 
of the constraints. However, if we are interested in accuracy as well as precision 
then any systematic differences between constraints or data should also be taken 
into account. In the most extreme case they could be in conflict to the extent that 
each of them would require its own set of parameter values, since they are in different 
regions of parameter space. This scenario can be considered as another hypothesis or 
model, i^O; to be compared with the other. Hi. Bayesian inference provides a very 
easy method of distinguishing between the scenarios Hq and Hi. This technique 
was illustrated in Refs. p^ , p| by examining the compatibility between B — )■ X^j 



branching ratios and the anomalous magnetic moment of the muon in the CMSSM. 

The natural logarithm of the ratio of posterior model probabilities provides a useful 
guide to what constitutes a significant difference between two models: 



log AE = log 



P{Hi\B) 
P{Ho\T>) 



log 



Zi P{Hi] 
ZoPiH 



0) 



(2.9) 



In Tab. |1| we summarise the conventions we use in this paper to separate between 
different levels of evidence. 



1 logAEl 


Odds, Zi/Zo 


Probability 


Remark 


< 1.0 


<3:1 


< 0.750 


Inconclusive 


1.0 


~ 3 : 1 


0.750 


Weak Evidence 


2.5 


~ 12 : 1 


0.923 


Moderate Evidence 


5.0 


~ 150 : 1 


0.993 


Strong Evidence 



Table 1: The scale we use for the interpretation of model probabilities. Here the log 
represents the natural logarithm. 
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The evaluation of the multi-dimensional integral (|2.2|) is a challenging numerical 
task. Standard techniques like thermodynamic integration |^2|] are extremely com- 
putationally expensive which makes evidence evaluation typically at least an order of 
magnitude more costly than parameter estimation. Some fast approximate methods 
have been used for evidence evaluation, such as treating the posterior as a multi- 
variate Gaussian centred at its peak (see e.g. but this approximation is clearly 
a poor one for multi-modal posteriors (except perhaps if one performs a separate 
Gaussian approximation at each mode). The Savage-Dickey density ratio has also 
been proposed |^ as an exact, and potentially faster, means of evaluating evidences, 
but is restricted to the special case of nested hypotheses and a separable prior on 
the model parameters. Bridge sampling allows the evaluation of the ratio 
of Bayesian evidence of two models and was implemented in fits to the CMSSM in 
Ref. [0. Bank sampling (a particular implementation of the MCMC technique) 
also allows evidence ratios to be calculated with half the number of points than re- 
quired for bridge sampling, and was used in CMSSM fits in Ref. []T3|. It is not yet 
clear how accurately bank sampling can calculate these evidence ratios, and both 
bridge and bank sampling have the disadvantages of only being able to calculate 
evidence ratios, not the evidence itself. When comparing models, this introduces 
an inefficiency of a factor of 2(A^ — 1)/N compared to the latter case. Various alter- 
native information criteria for model selection are discussed by [Q, but the evidence 
remains the preferred method. 



The nested sampling approach, introduced by Skilling [pi 
method targeted at the efficient calculation of the evidence 
terior inferences as a by-product. 



Feroz & Hobson |E9 



is a Monte Carlo 
It also produces pos- 
Ol built on this nested 



sampling framework and have recently introduced the MultiNest algorithm which 
is efficient in sampling from multi-modal posteriors exhibiting curving degeneracies, 
producing posterior samples and calculating the evidence value and its uncertainty. 
This technique has greatly reduced the computational cost of model selection and 
the exploration of highly degenerate multi-modal posterior PDFs. We employ this 
technique in this paper. 



3. The pMSSM 

The MSSM Lagrangian has the form £ = ^susY + '^soft where ^susy contains all of 
the kinetic terms, gauge and Yukawa interactions while preserving SUSY invariance. 
It is based on the gauge group G = SU{3)c x SU{2)l x U{\)y and superpotential, 
W , constructed with a particle content in the following chiral superfields shown with 
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their corresponding G charges: 



L: (1,2, -i), E: (1,1,1), 



Q:(3,2,i), f/:(3,l,|), 

D 6 



D : (3,1,-i), H^: (1,2,-^), H2 : (l,2,i). 
The superpotential is given by 



(3.1; 



(3.2) 



Here we use the convention in |^ and denote the SU{?>) colour index of the fun- 
damental representation hy x,y = 1,2,3; the SU{2)l fundamental representation 
indices by a, 6 = 1,2 and the generation indices by i,j = 1,2,3. eab = e"* is the 
totally antisymmetric tensor, with eu = 1. 

The soft SUSY-breaking part of the Lagrangian consists of different mass and 
coupling terms: 



soft 



-"gauginos I '^sfcrmions ~l~ '^trilinear ~l~ ''-'higgs 

where the part including the SUSY breaking sfermion masses is 



(3.3) 



■'sfermion 



,2n ~* 



(3.4) 

Each mass parameter in Eq. is a 3 x 3 hermitian matrix in generation space. 

-^higgs = mlH^^'Ht + ml^H2a*m + ^ab{mlH^Hl + H.c.) (3.5) 

gives the SUSY breaking higgs masses and bi-linear coupling terms. The SUSY 
breaking scalar trilinear couplings are 



"^trilinear ^ab 



H.c. 



(3.6) 

where fields with a tilde represent the scalar components of the corresponding capital 
letter superfield and the soft SUSY-breaking A-terms, each a complex 3x3 matrix 
in generation space, are defined (no summation on j is inferred) as 

{^u,D,E)ij = {Ua, Da, EA)ij/{Yu^D,E)ij- (3.7) 

Finally, writing the bino as b, w^='^'^^^ as the unbroken-S'?7(2)i gauginos and g^=^^-^^ 
as the gluinos, then the gaugino-mass part of the Lagrangian is 

3 8 



gaugino 



Mill + M2 ^ w^w"^ + M3 ^ g^g^ + H. 



c. 



A=l 



X=l 



(3. 



- 11 - 



The parameters together make a total of 105 free parameters in Lsoft, before re- 



phasing and higgs potential minimisation q^, In the CMSSM, the SUSY 
breaking scalar masses the gaugino masses and trilinear couplings are collapsed to the 
flavour independent parameters mo, M1/2 and Aq respectively, at grand unification 
scales Mgut ~ 2 x 10^^ GeV. m\ and are related to the Z— boson mass mz through 
higgs potential minimisation conditions. sign{^) = ± and tan /3, the ratio of the 
Higgs vacuum expectation values (vevs) remain as free parameters. However in this 
paper we instead explore the parameters at the weak-scale in its phenomenologically 
most relevant directions and following Ref. pTl call the set-up the pMSSM. In 



Subsection 3.1 we describes how the 25 parameters in the pMSSM set-up are derived 



from the much larger parameter space of the parent MSSM. 
3.1 Parameters 

Sources of CP-violation in the MSSM are tightly constrained by experimental lim- 
its on the electron and neutron electric dipole moments and from results on K-meson 
system experiments. Assuming that the MSSM soft SUSY-breaking parameters are 
real is consistent with such tight bounds, indeed significant departures from this 
assumption usually require a specific cancellation or suppression mechanism in or- 
der to pass the constraints. To suppress flavour changing neutral current (FCNC) 
processes, all off-diagonal elements in the sfermions masses and trilinear couplings 
are set to zero and the first and second generation soft terms are set to be equal. 
At, Af,, At- may all change the likelihood significantly, and we also include A^ = A^ 
because it is relevant for the computation of the anomalous magnetic moment of the 
muon We set Au = Ac = Ad = As = Q since these are proportional to the 

SM Yukawa couplings which are very tiny and so they will have negligible effect on 
the likelihood. All the other trilinear couplings are set to zero. Our Higgs-sector 
parameters are specified by w^l/J, and as discussed above, we must add tan/3 

and sign{fi) to the list of parameters. 

All of the parameters mentioned so far are purely non-SM. However, some of 
the SM parameters significantly affect the likelihood. The relevant SM parameters 
include the electromagnetic coupling constant aemi'm'z)^^^ and the strong coupling 
constant as{mz)'^^^ ■ The values of these two couplings are taken at the Z-boson 
pole mass rriz energy scale evaluated in the MS renormalisation scheme. The tau 
lepton mass and Gp, the Fermi constant have been so precisely determined that their 
uncertainty has negligible error on the likelihood and so they are fixed at their global 
average values: = 1.777 GeV and Gp = 1.16637 x 10~^ GeV~^p3|. The top and 
bottom quark masses are not as precisely known and can have significant effects on 
predictions of supersymmetric models. They are therefore included as parameters 
with, using the experimental measurements of their central values and uncertainties, 
Gaussian priors. Despite the fact that the Z-boson mass, mz, is precisely determined. 
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we include its uncertainty because one of the observables used for the analysis (the 
total decay width of the Z-boson Tz) is proportional to m| and the pMSSM predicted 
values can fall outside the expected experimentally determined range with a sigma 
variation in mz- So adding the SM parameters, 

OsM = {mz, mt, mbirrib)'^, aem{mzf^ , a,(mz)*^} (3.9) 

makes a total of 25 continuously varying parameters in the pMSSM. These are listed 
together with their ranges or Gaussian prior distributions in Tab. ^ These pMSSM 
directions make up its parameter-space. Our aim is to eventually construct a detailed 
map of the parameter-space that could be of help to or guide for collider and other 
SUSY-related experiments. In Subsection we briefly describe the observables 
considered and summarise the experimental constraints coming from each. 

3.2 Observables and experimental constraints 

The SM fits high precision electroweak data well However, on the one hand. 



there are some observables whose SM predicted values significantly differ from the 
corresponding experimental indications. The discrepancies could be explained by the 
direct or indirect presence of supersymmetric particles (or sparticles) in the inter- 
actions. On the other hand the very precise agreement between the SM prediction 
and the experimentally determined values of some other set of observables could be 
altered by the presence of non-SM particles. The absence of any significant such 
deviation puts tight constraints on possible new physics beyond the SM (SUSY in 



our case); see for instance |5^ and references therein. The values of the spar- 
ticle masses affects these tendencies. For instance, the effect of sparticles in loop 
corrections to electroweak physics observables (EWPO) decouple if their masses are 
much heavier than niz (300 GeV and above according to Ref. |^). Lighter spar- 
ticles with masses just above current experimental limits will alter the agreement 
between the electroweak data and SM predictions - hence the preference for low 
energy (weak-scale) SUSY. 

For the pMSSM set-up we use observables and constraints from high precision 
collider measurements that are sensitive to effects of new physics via virtual loops. 
These include five EWPO: the VT-boson mass, miy, the effective leptonic weak mix- 
ing angle, sin^^^y^, the total Z-boson decay width, Tz-, the anomalous magnetic 
moment of the muon, {g — 2)^ and the mass of the lightest MSSM Higgs boson, 
rrih] five B-physics observables: branching ratios BR{B — )• Xg'j), BR{Bs — )■ /x+zi"), 
BR{Bu- — )■ T~v), BR{Bu — K*'y) and the Bg mass-mixing parameter AM^j and 
the cosmological observable, dark matter relic density from WMAP5 results. We 
next briefly describe each of these physical observables and state the correspond- 
ing experimental constraints. We first discuss constraints from EWPO and end the 
Section with discussion of sparticle mass limits. 
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Parameter 


Description 






Prior range 


Ml 


Bino mass 






[-4 TeV, 4 TeV] 


Mo 


Wino mass 






[-4 TeV, 4 TeV] 


M3 


Gluino mass 






[-4 TeV, 4 TeV] 




1st /2nd generation Li slepton 


[100 GeV, 4 TeV] 




3rd generation L 


L slepton 


[100 GeV, 4 TeV] 




1st /2nd generation Er slepton 


[100 GeV, 4 TeV] 




3rd generation 


slepton 


[100 GeV, 4 TeV] 


777,- = 772,7 = 777- = 777,- 
Cl]^ *L 


1st /2nd generation Ql squark 


[100 GeV, 4 TeV] 


777r = 777r 


3rd generation Ql 


squark 


[100 GeV, 4 TeV] 


777- = 777r„ 


lst/2nd generation Ur squark 


[100 GeV, 4 TeV] 




3rd generation Ur 


squark 


[100 GeV, 4 TeV] 


777 7 = 777;„ 


lst/2nd generation Dr squark 


[100 GeV, 4 TeV] 




3rd generation Dr 


squark 


[100 GeV, 4 TeV] 


At 


Trilinear couplin 


g for top quark 


[-8 TeV, 8 TeV] 




Trilinear couplin 


g for b-quark 


[-8 TeV, 8 TeV] 




Trilinear couplin 


g for r-quark 


[-8 TeV, 8 TeV] 




Trilinear couplin 


g for /i-quark 


[-8 TeV, 8 TeV] 




up-type Higgs doublet mass 


[100 GeV, 4 TeV] 




down-type Higgs 


doublet mass 


100 GeV, 4 TeV 

L "J 


tan /3 


scalar doublets vevs ratio 


[2, 60] 




top quark mass 


54 


1 


172.6 ± 1.4 




Z-boson mass 






91.187 ± 0.021 




b-quark mass 






4.20 ± 0.07 




electromagnetic coupling constant 


127.918 ± 0.018 




strong coupling constant 


0.117 ± 0.002 



Table 2: The 25 parameters of the pMSSM model. The first twenty non-SM parameters are 
listed with their corresponding prior range. Gaussian priors are used for the SM parameters, 
which are the last five listed, along with their central values and standard deviations. 



VT-boson mass, mw' 

The CDF Run II electroweak public results cited the W^-boson mass measure- 
ment as the single most precise measurement to date and quotes [^| 

mw = (80.399 ± 0.025) GeV. (3.10) 

Theoretically the mass can be calculated from 

mw = — ;= 1 (3.11) 

^GHl-mym|)(l-Ar)' 

where aem is the fine structure constant at the mz renormalisation energy 
scale, Gp is the Fermi weak coupling constant and Ar includes all radiative 
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corrections to the mass (see e.g. 31 1 and references therein). The high 



precision in this measured quantity constrain any radiative corrections from 
new physics effects. The experimental precision is very high to the extent that 
measurements can be sensitive even to two-loop effects involving sparticles. We 
include a theoretical uncertainty of 10 MeV on mw by adding it in quadrature 
with the experimental uncertainty. We use SUSYPOPE ^ to calculate the 
VT-boson mass myy and the other EWPO. The most complete available SM 
two loop corrections and the dominant results for two loop SUSY corrections as 
implemented in SUSYPOPE currently give the most accurate predictions within 
the MSSM. 

Z-boson decay width, Tz'- 

The partial Z-boson decay width in the massless fermion limit {m'j/m^^ —¥ 0) 
is given by ||6^ 

rz^ff- = Nf ^ 6qcd{£^' + + ^e^/QCD (3.12) 

where g^^a are the neutral weak coupling constants modified to include elec- 
troweak (EW) radiative effects, 6qcd parametrises the QCD corrections and 
^ew/QCD includes some non-factorisable EW contributions. The colour factor 
NJ. is 1 for leptons and 3 for quarks. The current experimental value for the 



total decay width is 59 



Tz = (2.4952 ± 0.0023) GeV. (3.13) 

Theoretically 

rz = ri + rh + r^,. (3.14) 

where are the decay widths into SM leptons and quarks. Tinv is for the 
decays into invisible particles (neutrinos and possibly, if they are light enough, 
neutralinos) . Supersymmetric contributions enter via virtual corrections to the 
partial decay widths into lepton and quarks. 

Effective mixing angle, sin^6'gj'j: 

The effective electro-weak mixing angle depends only on the ratio of the effec- 
tive weak couplings 

Re{g,/ga) = l-Asm^eljf (3.15) 

for the vertex that couples the Z-boson and leptons / in the Lagrangian: 
ilj'^ [g^ — galb) Z^l. It is determined from various asymmetry measurements 
around the Z-boson peak from e+ e~ colliders after removing QED effects 
We use the experimental estimate |^ 



sin^ eljf = 0.2324 ± 0.0012. (3.16) 
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Z-pole asymmetry parameters from e+e — )■ // processes: 

The results from the LEP and SLC e~^e~ colhders on Z-boson properties (its 
mass, partial and total widths, and couplings to fermion pairs) are in good 
agreement with the SM predictions ||6^. The precision is high enough to probe 



loop-level predictions where both SM and beyond the SM corrections are ab- 
sorbed into effective coupling constants. The LEP data consist of hadronic and 
leptonic cross sections, leptonic forward-backward asymmetries, r polarisation 
asymmetries, bb and cc partial widths and forward-backward asymmetries. The 
Z-boson parameters derived from the data which we employ for our analysis 
include the ratios Ri (which we assume to be the average of Re, R^ and Rr), 
Rb and Rc. These are defined as 

_ TjZ ^ bb) _ TjZ ^ cc) T{Z^ m-) 

-Kh — „ , „ ; : r, JrCr — „ , „ ; ; r, it; 



T{Z hadronsY T{Z hadrons)^ T{Z ^ hadrons) 

(3.17) 

and are constrained to be 

Rb = 0.21629 ± 0.00066, R^ = 0.1721 ± 0.0030 and Ri = 20.767 ± 0.025 



The Z-boson interacts with fermions through a mixture of vector and axial- 
vector couplings. This makes the strength of the interaction between left- 
and right-handed fermions unequal and leads to the production of polarised 
Z-bosons at the e"'"e~ colliders. As result there are measurable asymmetries 
(such as a forward-backward asymmetry) in the angular distributions of the 
final-state fermions //. The forward-backward asymmetry is related to the 
probability that the / travels in the same (forward) or opposite (backward) 
direction to the incident e~ direction and is quantified by 

A,s = (3.18) 

where (JF{aB) is the cross section in the forward (backward) directions. In terms 
of the effective vector and axial- vector neutral current couplings, gvf and gAf 
respectively, other Z-pole asymmetries are: 

A^'f --A'Af af = l3vfgAi_ , . 

yvf 

Here A-^ gives a measure of the asymmetry for the different possible final state 
fermions. At LEP the Z-pole forward-backward asymmetries v4^^ and v4^^ 
were precisely measured for the final states bb and cc respectively. We impose 
A^^ = 0.0992 ± 0.0016 and A% = 0.0707 ± 0.0035 for our analysis. 
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Using polarised beams, the SLD experiment made a direct and precise mea- 
surement of the parameter Ae from the left-right asymmetry 



A,n = (3.20) 

were ai and are the e+ e~ production cross sections for Z bosons produced 
with left- and right-handed electrons respectively. The same parameter, A^, is 
also indirectly constrained by LEP experiments. Using the measurements of 
the parameters A^^ A'^ A^ and A"^ can then be inferred from ApB measurements 
at LEP. Hence, the LEP and SLC results form a complete set of the A-l" 
parameter measurements. The asymmetry parameter constraints we use are 
= 0.923 ± 0.020, = 0.670 ± 0.027 and A' = 0.1513 ± 0.0021 = A^. 

Muon anomalous magnetic moment, (5a^: 

The world average for the muon anomalous magnetic moment as determined 
from e'^e~ — )■ hadrons-based experiment at Brookhaven |^6|] is a^^^ = ^{g — 



2)^ = 1165920.80 ± 0.63 x 10 . Other results from experiments based on 



the r lepton decay to hadrons [0, ^ ^ differ slightly and we are not using 
those here^. The experimental results are around 3cr deviation from the SM 
prediction |7T[ aj^^ = 1165917.85 ± 0.61 x 10"^ giving 



Sa^ = a'^^P - a^*^ = 29.5 ± 8.8 x 10"^°. (3.21) 

SUSY is a good new physics candidate that can explain this deviation. For 
sparticles all of mass Msusy, their contribution to is of order 

af ^ 13 X 10-1° ( tan/3sign(^). (3.22) 
V Msusy J 



We use micrOMEGAs |]73|, 0, |75|, |76| to predict a^, correcting it with A^, which 
includes contributions from the dominant QED-logarithmic piece, some domi- 
nant two-loop corrections and the recently computed tan /9-enhanced term 
(see 



77|| and references therein) 



^susY ^ ^susY,iL f 1 _ 1^ log M^] ( _J_) + 2-loop terms. (3.23) 

For our analysis we add in quadrature a theoretical error of 2.0 x lO^^*^ to the 
8.8 X 10"^'^ error above and use the average value (5a^ = (30.2 ± 9.2) x 10"^". 

Lightest Higgs boson mass, rrih: 

The SM higgs mass is constrained to rrih > 114.4 GeV by LEP direct search 
experiment (see e.g. [|7S| and references therein). The predicted higgs mass can 



^See IfTOi for a recent review. 
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be parametrised by ghzz / dhzz^ ratio of the MSSM higgs coupling to two 
neutral Z-bosons to the equivalent SM coupling. In the MSSM gnzz/dhzz — 
sin^(/3 — a) so we use the above LEP mass limit for sin^(/3 — a) > 0.95 and use 
the MSSM higgs mass limit nih > 89.7 for sm^{/3 - a) < 0.95. 

Branching ratio (B Xg'j): 

The experimental value of the decay rate of the flavour changing process 
B — )■ Xg'y agrees to high precision with the SM prediction. This stringently 
constrains new physics models that may contribute to the process. For SUSY 
the leading contributions come from loops with charged Higgs bosons (these 
interfere constructively with the SM contributions) and charginos. Loops with 



neutralinos are small (see for instance, 17911). A recent theoretical estimate 



for the SM contribution to the branching ratio of S — )■ which we call 



BR{b S7) at next-to-next-to leading order (NNLO) in QCD is |80|, gl], |2|, |3 
BR{b s-f)sM = (3.28 ± 0.23) x 10"^, where a 1.6 GeV lower energy cut is 
applied to the photon. The central value is ~ la lower than a world-average 
experimental value from the Heavy Flavor Averaging Group ||8^ (HFAG)^ 



BR{B Xs-f)exp = (3.55 ± 0.22V^fQ ± 0.03) x 10~^ (3.24) 

The values predicted for this observable at different parameter points are con- 
strained by combining in quadrature the experimental and SM prediction er- 
rors: 

BR{B X,j)exp = (3.55 ± 0.42) x 10"^ (3.25) 

The branching ratio therefore constrains non-SM contributions. We use Su- 
perlso2.0 |^ to predict the MSSM plus SM branching ratio®. 



Branching ratio {Bs — t- yU^/x ): 

The branching fraction for the flavour changing process Bs is predicted 

to be (3.42±0.54) x 10"^ within the SM In the MSSM, interactions involv- 



ing neutral Higgs bosons can enhance the branching fraction by several orders 
of magnitude at high tan /3. The branching ratio is experimentally bounded 
from above by recent CDF II results, implying 

BR{Bs fJ.^fi') < 5.8 X 10"*^ (3.26) 



at 95% CL [^. We apply a continuous likelihood constraint derived from 
CDF II data upon the MSSM prediction ||89[. The resulting penalty is shown 
in Fig. |I|. 



^We note that the most recent central value has shifted slightly to 3.52: this small change would 
leave no significant imprint on our fits. 

^The program SusyBSG [^6| now contains a more accurate prediction of this branching ratio to 
two-loops in MSSM parameters. It could lead to an estimated shift in the predicted branching ratio 
of 0.13 X 10^"*, which may be considered to be included within our estimate of theoretical error. 
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-8 -7.8 -7.6 -7.4 -7.2 -7 
logio[BR(B^->^V)] 



Figure 1: Likelihood penalty on the predicted value of BR{Bs — )■ /^^/U ). 



Bg-Bg mass difference, AMb/- 

The neutral i?-meson oscillates between particle and antiparticle states via 
flavour changing processes. The frequency of oscillation is proportional to 
the mass difference AMb^ It has been measured to be {AMbJexp = 17.77 ± 
0.12 ps~^ [Q. Its SM prediction can be obtained from an overall unitarity tri- 
angle fit: {AMbJsm = 20.9±2.6ps~^ ||9l|. We use the ratio of the experimental 
constraint to the SM prediction 

to constrain the predicted frequencies. Our pMSSM predicted values are based 
on the results in 11921 and references therein: 



TjMSSM {.^Mb,)mSSM ,o oo\ 

[IXMbjsm 

^ 647r sin^ 6^ mi,{^i,)ms{^b) (ey tan^/3)^ 

ae^MlSoimj/ml,) [1 + (eo + ey)tan/3]2[l + eotan/3]2 



where m^^s are the bottom- and strange-quark masses evaluated at the MS 
scale fj,b = 



are the effective couplings that parametrise the correction to the down-type 
Yukawa couplings, 

uf \ xlogx yhgy 
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/i is the supersymmetric Higgs mass terms and At is the trihnear soft breaking 
term involving the stops. 5*0 is given by 



Sf)(x) 



Ax — llx + x^ 
4fl - x)2 



3x^ log X 
2(1 -x¥ 



Branching ratio {B^ tu): 

The purely leptonic decay B^- — )■ r~z/ proceeds via the annihilation of b- and 
"U-quarks into W~ . The SM prediction for the branching ratio of the process is 
given by 



BR{B^ Tv) 



SM 



57r 



ml 



m 



B J 



(3.31) 



where tub and are the B meson and r pole masses, respectively, and tb 
is the i?~-meson lifetime. For the SM prediction we use the average of the 
result from unitarity triangle fits {BR{Bu — ?■ tu) = 0.85 ±0.14 x 10~^) and 
the result obtained from the experimental determination of Vub and /b^/Bb^ 
{BR{Bu — > Tu) = 1.39 ± 0.44 x 10~^) adding the errors in quadrature to: 



SM 



1.12 ±0.46 X 10" 



(3.32) 



. For the experimental constraint upon the branching ratio, we use the average 
of the Belle and BaBar experiments, adding their errors in quadrature: 



BR{Bu rz/)e^p = 1.41 ± 0.43 x 10" 
Eqs. ( |3.32D , ( |3.33D are then used to form the constraint 

BR{Bu Tu)exp 



1.26 ±0.41. 



BR{Bu Tiy)sM 
For the pMSSM predictions we follow |^ and predict 



(3.33) 



(3.34) 



R 



MSSM 
Btv 



BR{B^ Tv] 



MSSM 



BR{Bu rv] 



SM 



m 



tan^ (3 



mjj± / 1 ± eotan/3 



(3.35) 



where eo is the effective coupling defined in Eq. ^.29| , m^^ is the the B-meson 
mass and mj^± the charged higgs boson mass. 

Aq-: B — )■ K*'y Isospin asymmetry: Isospin symmetry predicts the ampli- 
tudes for the decays B^ — )■ K*^'y and B~ — )■ K*~'y to be equal at leading order 



in perturbation theory. Isospin-breaking effects in the process B — )■ K*'y [p3 | 



^See "New Constraints from B and K rare decays" at http://www.utfit.org/ 
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may therefore provide a sensitive probe of physics beyond the SM. The isospin 
asymmetry for the exclusive process B — )■ K*'y is defined as 

. _ TjB^ ^ K*^^) - r{B- K*'^) 



The world average experimental value is |£3| 

Ao- = 0.0375 ± 0.0289. (3.37) 



In order to fit this, we use the MSSM prediction from Superlso2.0 [^] which 
contains the available next-to-leading-order (NLO) contributions to Aq-, in- 
cluding the complete supersymmetric QCD corrections to Wilson coefficient 
operators. It also includes some partial next-to-NLO (NNLO) SM QCD cor- 
rections. 

• Dark Matter relic density: 

The Wilkinson Microwave Anisotropy Probe (WMAP) fits to a cosmological 
constant plus cold dark matter model (ACDM) imply a dark matter relic den- 
sity of VLoMh^ = 0.1143 ± 0.0034, where h is the reduced Hubble constant 
P^ . We assume i?-parity, resulting in a stable lightest supersymmetric parti- 
cle (LSP). The neutralino, Xi? LSP then has the correct properties to make up 
the cold dark matter, being massive, stable and neutral and we constrain the 
prediction of its relic density to lie on the WMAP5 central value, but inflate 
the uncertainties with an assumed error on the theoretical prediction: 

nDMh"^ = 0.1143 ± 0.02. (3.38) 

We summarise the experimental constraints used in our fits in Tab. ^, listing the 
relevant references with each. 

3.3 Direct search mass limits 

The absence of sparticle or Higgs boson production at current collider searches 



for supersymmetric particles puts lower bounds on their possible masses We 
veto any pMSSM points that violate the limits. The limits are derived from various 
experiments that usually a priori assume the validity of a chosen model (usually 
the CMSSM). Where possible, we use more appropriate model independent limits 
upon sparticle masses coming from searches. SUSY particles may be pair produced at 
colliders that have sufficient energy, then undergo subsequent decay into SM particles 
and neutralino LSP. Hard jets or leptons associated with missing energy coming 
from the neutralino then constitute SUSY direct search signatures. Constraints on 
sparticle of mass m often dependent upon the mass difference Am = m — rnisp 



which is correlated with the energy of visible sparticle decay products | 106 , [2^ . The 
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Observable Constraint 



Th. Source 



Ex. Source 



mw [GeV] 
Tz [GeV] 

Qlep 



sm 



eff 



6a^ 

Rl 
Rl 

^FB 
^FB 

= A" 

A' 

A^ 

Br{B ^ X,7) 
Br{Bs /i^/i"^ 
Ramb^ 

RBr{Bu~^Tu) 



80.399 ±0.027 
2.4952 ± 0.0025 
0.2324 ±0.0012 
(30.2 ±9.0) X 10^° 
20.767 ±0.025 
0.21629 ±0.00066 
0.1721 ±0.0030 
0.0992 ±0.0016 
0.0707 ±0.035 
0.1513 ±0.0021 
0.923 ±0.020 
0.670 ±0.027 
(3.55 ±0.42) X 10^ 
see Fig. |1| 
0.85 ±0.11 
1.26 ±0.41 
0.0375 ± 0.0289 
0.11 ±0.02 



m 



|63 



[96 



98, 99 



m 



m 



|62 



m 



m 



m 



80, BTl B2, 100 



|73|, ^0^1 



m 



m 



95 



m 



35 



35| 



3^1 



35 



35 



35 



B8 



90 



m M ITU3I 



m [HI 

^, ^ ^ [HI 



Table 3: Summary of constraints on predictions. Theoretical uncertainties have been added 
in quadrature to the experimental uncertainties quoted. 



sparticle mass limits derived then depend on this energy; depending upon whether 
Am is low (5 - 10 GeV) or higher. The mass limits we impose on the sparticles and 
the lightest Higgs boson are summarised in Tab. ^ |]S3 . 



A recent random scan study of the pMSSM [Q found that a CDF/DO bound ||107 



108|| was quite restrictive on their pMSSM parameter-space random sample points 
when the relic density constraint was only applied as an upper bound (i.e. allowing 
for additional extra-MSSM sources of dark matter). The bound states that, for 
m^± - m^o < 50 MeV, m-± > 206 If/i^P + 171 |f/i/,P GeV at 95% confidence level 

Al Al Xl ' ' ' ' 

. Here \Uih\ are the Wino and Higgsino content of the lightest chargino, 

respectively. We did not use this bound in our MultiNest sampling procedure, but 
we have checked retrospectively that it would not have significantly changed our fits 
since only less than 1% of the posterior probability density fails this constraint. 



3.4 Sampling procedure 

In this Subsection we summarise the Bayesian inference elements (briefly reviewed 
in Section ^ for the pMSSM and the sampling procedure we employ for fitting it to 
the indirect collider and cosmological data. All of the pMSSM parameters, 6, listed 
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condition sparticle mass lower bound/GeV 



sin2(a - /3) > 0.95 mt 114 

sin2(a - (3) < 0.95 rrih 89.7 

rrifj — m^o > 10 GeV m^^ 87 

rrif-^ — m^o < 10 GeV m^^ 73 

mg^ — m^o > 10 GeV m^^ 100 

mg^ — m^o < 10 GeV mg^ 73 

m^^ — m^o > 10 GeV m^^ 95 

m^^ — m^o < 10 GeV m^^ 73 

'^i^^ — > 10 GeV 94 

— m^o < 10 GeV mj>„ 43 
m;>^ > 300 GeV m-± 43 

< 300 GeV m-1 92.4 

m^o 50 

— m^o > 10 GeV 95 
mj — m^o < 10 GeV 65 

— m^i) > 10 GeV 95 
rrij — m^o < 10 GeV mi 59 

mg 318 

mg 195 



Table 4: T/ie lower bounds applied to MS SM particle masses. 



in Tab. H, are varied simultaneously, our calculation being driven by MultiNest and 
the high energy physics software mentioned in the following paragraphs. MultiNest 
is described in Appendix 0. Each point is passed in SUSY Les Houches Accord 
(SLHA) format ||109|| to the different particle physics software we use to predict the 



observables described in Subsection For each set of parameters MultiNest 
selects, the following steps are followed: 



1. The (input) parameters are passed to S0FTSUSY2 . . 18 |^ which produces 



the MSSM sparticle masses and couplings. Unphysical points are flagged by 
the program to have one or some combination of the following properties: ab- 
sence of electroweak symmetry breaking, the presence of a tachyon, a non- 
perturbative point where the calculation can no longer be trusted or the lack 
of numerical convergence (which usually occurs close to a boundary of good 
electroweak symmetry breaking) . If any of these properties are flagged then the 
point is discarded before any further computations and the parameter point is 
given a zero likelihood. In addition to this, sparticle spectra that violate mass 
limits, shown in Tab. |[ from sparticle searches or that have a non-neutralino 
LSP are also discarded. 
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2. Physical parameter points are passed in SLHA format to micrOMEGAsS . 2 
[f^ , [75| , |76| which calculates the neutralino dark matter relic density, the branch- 
ing ratio BR{Bs — )■ /^"""/x") and the anomalous magnetic moment of the muon 



5a^. The physical point is then passed to the computer program SuperIso2 . |^ 
and other subroutines. The former calculates BR{B — )■ Xs'~f) and the isospin 
asymmetry in B meson decays Aq_ while the latter computes the B-physics 
ratios Rbu^tu and Ramb ■ 



3. We then use SUSYPOPE |^ ^ to calculate the VT-boson mass mw, the effective 
leptonic mixing angle variable sin^^^^^, the total Z-boson decay width Tz and 
the other electroweak physics observables listed in Tab. ^ to two loops. 



The various physical observables described in Subsection p.2| , which are derived 
from the input parameters, form the data set, Di, i = 1, ... ,19: 

D = {mw, sin' e[j^, Tz, 6a„ R^ A%, A' = A\ R^, A%, A"''^, (3.39) 
BR{B ^ Xs-i), BR{Bs ^ ii'), Ao-, Rbr(Bu^tu), Ramb,, 

For each predicted value Oi of observable Di, the corresponding likelihood is 



p(A|e,iJ) = -i=exp 



where /Xj and cxj are the experimental central values and errors given in Tab. ^. We 
assume that the observables are independent and combined the likelihoods to 

19 

L(0) = P(D|0, H) = \{ P(A|e, H) LBRiB.^,^,-) (3.41) 

1=1 

where LBR{Bs^^l+^l-) is the likelihood value for the indicated observable (see Fig. |l]). 
The predictions from the physical points, as enumerated above, are checked against 
experimental values and the deviations from this are quantified by the individual 
likelihood functions Eq. p.40| . The likelihoods from the different observables are 



combined into one overall likelihood, Eq. p.41| , which is then multiplied by the prior 



probability density Eq. to produce posterior probability density Eq. 

3.5 Computer resources 

We end this Section with the presentation of an estimate of the computing resources 
used. In the MultiNestI.3 algorithm we used 4,000 live points (see the Appendix 
and p9|, ^ for details) and more than 8.6 x 10^ likelihood evaluations for the linear 
prior and more than 2.1 x lO'^ for the log prior case. The overall number of likelihood 
evaluations in this work is more than 2.5 x 10^. The computing was performed by 
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79 12- hour jobs at the Darwin cluster of the high performance computing service 
(HPC) at the University of Cambridge. Each run used 128 threads or CPUs. We 
also used around 40 8-hour jobs at COSMOS (the UK Cosmology supercomputer at 
DAMTP). At COSMOS each run used 64 CPUs. The total run time adds to more 
than 15 CPU years. However the efficiency of the more recent version of MultiNest 
has improved and these computations could take about half of the stated time. 



4. Results and Analysis 

We now present the results of the global fit for the weak-scale pMSSM param- 
eters to indirect collider and cosmological data. We first give the marginalised 1- 
dimensional posterior probability distributions for the 25 parameters and sparticle 
masses in Section ^?1|. We also give, see Tab. ||, a good fit point in the pMSSM param- 
eter space (i.e. the sampled point with maximum likelihood) and the corresponding 



sparticle spectrum in Tab. [4.2| . The observables posterior PDFs are discussed in Sec- 
tion [4.3| . The values which each of the observables take, together with corresponding 
deviations from experimental values, at the good fit point are given in Fig.s ^ and ^ 
In Section [4.4| we present an implementation of Bayesian model selection by compar- 
ing two different pMSSM hypothesis each with either sign of the /i parameter. We 
applied the 'gaugino code' prescription to the pMSSM in Section ^75[ This is con- 
nected with the gluino-neutralino mass ratio and its role in discriminating different 
SUSY-breaking phenomenology models. We address the case of fine-tuning in the 
pMSSM parameters in Section I 



Results and analysis on neutralino dark matter and its direct detection prospects 
are presented in Section ^ A general feature of the results is that they exhibit some 
amount of prior dependence. This is clearly shown for the case of gluino-neutralino 



mass ratio in Subsection [4.5| , the amount of fine-tuning in the parameters in Subsec- 
tion [4.6| , the muon anomalous magnetic moment (see the plots in Fig. |]), and the 
dominant neutralino dark matter annihilation/co- annihilation channels addressed in 
in Section ^. It is expected that with more precise and direct data from the Teva- 
tron and/or the LHC would lift the prior dependence. However interestingly enough, 
some of the results, such as the mass of the lightest CP-even Higgs boson mass, are 
very similar for the different priors considered. 

4.1 Parameters and sparticle mass posterior PDFs 

Our assumptions about the pMSSM priors are presented in Section The 
priors are updated using the set of different experimental measurements explained 
in Section |3.2| via Bayes' theorem (see Section |2.1| ). A measure of the amount of 
information in the likelihood can be seen from the results of the sampling procedure 
in the form of posterior PDFs. One dimensional marginalised posterior PDFs of 
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the parameters are shown in Fig. ^. Most of the scalar mass terms in the log prior 




Figure 2: Marginalised 1- dimensional posterior PDFs for the pMSSM parameters. The 
soft-breaking scalar mass terms in the log prior scenario (broken lines) are mostly reduced 
with respect to the corresponding mass term in the linear prior scenario (solid lines). All 
masses are in TeV units and niq^ = m^^ = m^^ = = m^^. 



scenario are much reduced with respect to the corresponding mass term in the linear 
prior scenario, as expected. Changes to the scalar mass priors (in the log prior 
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case) do change the posterior PDFs of Mi and M2. This is because the dark matter 
hkehhood is driving the fit, and dark matter co-annihilation requires the mass of the 
hghtest neutrahno (controlled by the smaller of Mi and M2) to be close to that of 
the scalar that it is co-annihilating with (see Subsection p.2| for discussion on the 
co-annihilations in the pMSSM). On the other hand, M3 shows approximate prior 
independence. For linear priors, upper bounds on the parameters come only from our 
constraint upon the prior range, there not being sufficient power in the data yet to 
constrain the parameters from becoming very large. Mi, M2 and M3 have PDFs that 
tend to zero as the parameter tends to zero, but a finite bin-size means that this may 
not be evident in plots. Since ttihi and may take negative values, no logarithm 
was applied to them in the log priors case. Thus the apparent prior independence of 
rriHi is no surprise. The A terms show a peak structure because if their magnitudes 
are too large, the squark mass becomes tachyonic and thus disallowed. Radiative 
corrections to the lightest CP-even higgs mass involving the stops imply that they 
must be quite heavy, above ~ 2 TeV in order to push the higgs mass above its direct 
search bound. Thus, the log prior only disfavours points at m^^ = = 4 TeV by a 
factor of 1/2 compared to the lower values. In contrast, other sfermion masses which 
may be as low as 200 GeV obtain a log prior suppression of 1/20 at 4 TeV masses, 
i.e. the difference between log and linear priors is much more evident in that case. 
Thus there is no large enhancement of the posterior for small m^^ = mi^, which 
shows approximately identical posterior PDFs for the two different prior choices. 

We display the posterior PDFs of the pole sparticle masses in Fig. ^ along with the 
posterior PDF of the fi parameter; and the 68% and 95% Bayesian credibility ranges 
for the sparticle masses in Tab. ^. The mass distributions can be understood from 
the parameters PDFs discussed above, since there is a rough one to one correspon- 
dence between the tree-level mass and a mass parameter for many of the sparticles. 
For the third family sparticles, the A-terms and n parameter contribute via the large 
mixing, which is proportional to the analogous SM fermion mass. The first and sec- 
ond generation sfermion masses are approximately degenerate, since the degeneracy 
is only broken by terms of order the second generation fermion mass, negligible com- 
pared to the sfermion masses. We see approximate prior independence of m^^, m^^, 
rrig, m^^ and itia ~ tuh ~ m/^±. Other sfermion masses show the expected log prior 
suppression at high masses. 

The MSSM lightest CP-even higgs mass at tree level is, 

ruh = mz\cos2P\, (4.1) 

but it receives large radiative corrections from third generation particles, of order 
30% of its mass. Since the tan /3 posterior is similar for both priors, then the tree- 
level value of rrih will also be similar. An additional small prior dependence comes 
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Mass [TeV] 



pMSSM fit (linear) 
68% interval 95% interval 



pMSSM 
68% interval 



fit (log) 
95% interval 



m-g 



(0.120,0.124) 
(0.111,0.338) 
(2.192,4.000) 



0.117 
0.608 
1.280 



0.129) 
0.405) 
4.000) 



0.119 
1.129 
1.844 



0.123) 
3.284) 
4.000) 



(0.118 
(0.509 
(0.924 



0.128) 
3.996) 
4.000) 



m^o 
m^o 
m^o 
m^o 



(0.899,1.125) 
(0.772,1.687) 
(0.126,2.000) 
(0.165,3.571) 
(0.765, 1.428) 
(1.335,3.217) 



0.302 
0.322 
0.446 
0.649 
0.315 
0.648 



1.599) 
2.630) 
2.000) 
3.877) 
1.730) 
3.771) 



0.100 
0.301 
0.539 
1171, 
0.290 
1.174 



0.273) 
1.675) 
2.000) 
2.967) 
1.692) 
2.981) 



(0.100 
(0.176 
(0.237 
(0.585 
(0.162 
(0.586 



0.397) 
2.491) 
2.000) 
3.631) 
2.532) 
3.656) 



mr,. 



mi 



mi 



(1.901,4.000) 
(2.039,4.000) 
(2.124,4.000) 
(1.718,2377) 
(2686,4.000) 



0.664 
1.155 
1.230 
1.048 
1.832 



4.000) 
4.000) 
4.000) 
3.188) 
4.000) 



0.279 
1.000 
0.588 
1.000 
1.166 



0.398) 
0.606) 
1.392) 
0.884) 
4.000) 



(0.207 
(1.000 
(0.256 
(1.000 
(0.502 



0.922) 
2.196) 
3.149) 
1.965) 
4.000) 



ma 



m 



t2 



mr 



m 



b2 



(1.965,4.000) 
(1.634,2.841) 
(2.591,3.710) 
(1.432,2.679) 
(2.319,3.637) 



1.143 
1.233 
2.078 
0.952 
1.746 



4.000) 
3.423) 
3.968) 
3.320) 
3.967) 



0.590 
1.529 
2.424 
0.478 
1.891 



0.987) 
2.633) 
3.644) 
1.695) 
3.452) 



(0.434 
(1.105 
(1.912 
(0.347 
(1.354 



2.573) 
3.271) 
3.949) 
2.701) 

3913) 



Table 5: pMSSM sparticle mass ranges corresponding to 68% and 95% Bayesian credibil- 
ity, "^[^i^j represents mass of 1st and 2nd generation left-handed(right-handed) sleptons. 



mainly from m^-^ and m^^ prior dependence, but really the model itself constrains 
the higgs masses to be largely prior independent. The approximate mass degeneracy 
in the heavy Higgs masses and little dependence on priors can be seen to originate 
from the relationships between them. The tree-level MSSM pseudo-scalar CP odd 
higgs mass is given by 

+ <2' (4-2) 



2|/i|2 + m2 



and it receives quite large radiative corrections. There is only a small prior depen- 
dence of m^, and Eq. [4.2| shows that the approximate prior independence of m^ is 
something of an accident since m\j^ and yU show some prior dependence, but this 
largely cancels in its effect on mA- Notice that 



m 



H 



m. + my + 



m 



z) 



+ Am\m\ sin^(2/3) 



m 



mA-\-m 



w 



(4.3) 
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and since tha is usually far greater than mw and mz then tua rriH ^ f^H^ holds, 
although loop corrections contribute to a small non-degeneracy, which however is not 
visible to the eye. 

4.2 Good-fit points 

For a given prior, the point with the highest likelihood is termed the best-fit point. 
We estimate a 'good-fit point' by picking the point with the highest likelihood out 
of the many millions sampled by the MultiNest algorithm. We shall see below 
that our under constrained fits exhibit parameter degeneracies, and so we should 
not expect the particular values of the good-fit parameters to be unique. Monte 
Carlo sampling methods are also ill suited to finding best-fit points, which are much 
better determined by hill climbing algorithms, for example. On the other hand, the 
properties of such a point are instructive to view in order to get a feel for which 
observables are pulling in which direction. The values of the parameters at the good- 
fit points are given in Tab. ^ The good-fit point of the linear prior has quite large 
SUSY breaking mass parameters, of the order of a couple of TeV in most cases, 
and large tan/3. On the other hand, the log prior sampling found a good-fit point 
with a mixture of sub- TeV and several TeV masses, because it spent more time 
exploring points with lower sfermion masses. It also has a moderate value of tan/3. 
Interestingly, both good-fit points have a negative /i parameter. This is in contrast 
to the CMSSM case, where the anomalous magnetic moment of the muon prefers 



positive /i. We shall discuss the sign of /i further in Subsection below. 



We see the resulting MSSM spectra in Tab. [4.2|, where the log prior good-fit 



point has some light sleptons and squarks. For the linear prior, all of the sparticles 
have masses around the TeV scale or higher. 

The statistical pulls of the various observables are shown in Figs. ^ and ^ We see 
that the observable that most discriminates between our two good-fit points is the 
anomalous magnetic moment of the muon, which is much better fit for the log prior 
point, where it receives large corrections from lighter slepton and gaugino masses. 
171]^ and Tz also show a large difference between the two points, whereas most of the 
other observables display only a small difference in statistical pull between these two 
points. In fact, the log prior good-fit point has a significantly lower value (and 
therefore, a better likelihood). Were nested sampling perfect for finding the good 
fit point, the likelihood (or equivalently, the total x^) values would be the same for 
both priors. This difference is therefore due to the finite sampling resolution of the 
posterior distributions. 
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Parameter 


Description 


Linear prior nt 


Log prior tit 


Ml 


Brno mass 


-2947 


-250 


71 /f" 


Wmo mass 


-1297 


-3017 


71 /T 

M3 


Glumo mass 


-2397 


-642 




lst/2nd gen. Ll slepton 


1040 


174 




3rd gen. Ll slepton 


2640 


993 




lst/2nd gen. hji slepton 


2301 


201 




3rd gen. Er slepton 


3748 


3530 




lst/2nd gen. Ql squark 


878 


165 


mi^ = 


3rd gen. Ql squark 


2301 


0001 

2321 




Ist/ind gen. Ur squark 


3027 


1515 




3rd gen. Un squark 


2618 


2905 




lst/2nd gen. Dr squark 


1368 


329 




3rd gen. Dr squark 


iOo4 


126o 


At 


top quark trilinear 


1 CiC\0 

-iyod 


DDi 


A-b 


b-quark trilinear 


Q c; /1 1 
-o04i 


/2 / 


A 


7 -qudiK tiunieai 




oiyo 




//.-quark trilinear 


2154 


2951 




down-type Higgs doublet 


2548 


3445 




up-type Higgs doublet 


882 


669 


tan^ 


Higgs vevs ratio 


52.0 


21.0 


mt 


top quark mass 


173.4 


175.3 


mz 


Z-boson mass 


91.186 


91.190 




b-quark mass 


4.16 


4.26 




e- coupling constant 


127.95 


127.91 




s-coupling constant 


0.1168 


0.1161 


l^{Msusy) 


higgs parameter 


-942 


-770 



Table 6: pMSSM input parameters for both good-fit points. Mass parameters are in units 
of Ge V. The magnitude of the fi parameter is a derived quantity. 
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Figure 3: Marginalised 1-dimensional posterior PDFs for the fj, parameter and pMSSM 
sparticle masses, all in TeV units, for log priors (broken lines) and linear priors (solid 
lines). 1st and 2nd generation left-handed sleptons, left-handed squarks and right-handed 
squarks are collapsed into single parameters, mr , ruq^ and mr respectively. 
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Table 7: Tabular and pictorial representation of the sparticle masses at the pMSSM good 
fit point for both linear and log priors. All numbers in the table are in GeV units. The 
pictures on the right hand side shows the difference in the pattern of sparticle masses for 
the linear (upper picture) and the logarithmic (lower picture) prior case. 
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|Qmeas_Qf,t|^^„ 



Observable 



Measurement 



pMSSM fit (linear) 




[GeV] 
Tz [GeV] 

5a^x10'° 



A"=A' 

Ab 



A" 

BR(B^X3Y)x10'' 

u 



80.399 ± 0.027 80.370 

2.4952 ± 0.0025 2.4952 

0.2324 ± 0.0012 0.2314 

30.20 ±9.02 7.15 

20.767 ± 0.025 20.761 

0.21629 ±0.00066 0.21605 

0.1721 ±0.0030 0.1722 

0.1513 ± 0.0021 0.1480 

0.923 ±0.020 0.935 

0.670 ±0.027 0.668 

0.0992 ±0.0016 0.1038 

0.071 ± 0.035 0.074 

3.55 ± 0.42 3.63 

1.26 ±0.41 1.00 

0.85 ±0.11 0.99 

0.0375 ± 0.0289 0.0764 

0.11 ±0.02 0.10 




b h 2 B 



Figure 4: Statistical pull of various observables at the good-fit point from the linear prior 
sampling. The good-fit value and the measurement for each observable are listed. The 
statistical pull, defined to be the number of sigma the point differs from the experimental 
central value is shown by the horizontal bars. 
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|Qmeas_ Qfit| ^^meas 



Observable 



Measurement 



pMSSM fit (log) 



[GeV] 
Tz [GeV] 

5a^x10'° 



A"=A' 

Ab 



A" 

BR(B^X3Y)x10'' 

'^BR(B ^Tv) 

U 



80.399 ± 0.027 80.402 

2.4952 ± 0.0025 2.4964 

0.2324 ± 0.0012 0.2314 

30.20 ± 9.02 26.74 

20.767 ± 0.025 20.760 

0.21 629 ± 0.00066 0.21 962 

0.172110.0030 0.1723 

0.1513 ± 0.0021 0.1483 

0.923 ±0.020 0.935 

0.670 ±0.027 0.685 

0.0992 ±0.0016 0.1040 

0.071 ± 0.035 0.074 

3.55 ± 0.42 3.42 

1.26 ±0.41 1.00 

0.85 ±0.11 1.00 

0.0375 ± 0.0289 0.0748 

0.11 ±0.02 0.13 





2 B 



Figure 5: Statistical pull of various observables at the good-fit point from the log prior 
sampling. The good-fit value and the measurement for each observable are listed. The 
statistical pull, defined to be the number of sigma the point differs from the experimental 
central value is shown by the horizontal bars. 
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4.3 Observables' PDFs 



The posterior PDFs for the observables used to constrain the pMSSM are given 
in Fig. 1^. Differences between the two prior cases are mostly due to the fact that 
the sparticle mass PDFs are larger in the linear prior, leading to a suppression of 
SUSY effects in the loops of most observables. The -B-physics observables tend to 
have similar posterior PDFs in the two prior cases. In most of the EWPO, there are 
larger differences in the posteriors. For 5a^, there is a particularly large difference 
between the central values of the linear and log prior PDFs. The leading one-loop 
gaugino contribution at large tan {3 is given by ||9^ 



m^utan/3 , „ n x . . 

5a, ^ {glM.Fr + glM^F^) , (4.4) 

where Fi and F2 are positive loop functions proportional to m~^gy for the case of de- 
generate sparticles in the loops. The dominant contributions coming from gauginos 
and sleptons therefore lead to an enhanced value of Sa, when they are lighter, as is 



evident for the log prior fits. As we shall see in Section |5.2| and Tab. 0, the linear 
prior fits prefer Higgsino exchange to be the dominant LSP annihilation process, as 
opposed to slepton co-annihilation in the log prior case. This occurs at heavier neu- 
tralino LSP masses, and hence heavier smuon masses (which are always constrained 
to be heavier than the neutralino LSP). 6a is then relatively badly fit as can be seen 
in the good fit point example where the statistical pull is more than 2cr. 

4.4 Sign(/i) comparison 

The posterior PDFs in Fig. || indicates that the pMSSM prefers /i < compared to 
/i > 0. This is interesting since in the previous studies of CMSSM /i > was seen to 
be preferred by the combination of BR{b — 57) and (5a^. One of the statistical tests 
(a predictive likelihood ratio test) in Ref. found that the two measurements 
are incompatible in the CMSSM, but the other found no strong evidence for this 
and so the final conclusion of the analysis remains unclear. The sign of the SUSY 
contribution to BR{h — )■ 57) is dependent upon the sign of /i. There are two dominant 
SUSY contributions to consider: the first comes from diagrams involving a charged 
Higgs boson and up-type quarks. The second, involving a chargino and up-type 
squarks, depends upon the sign of the product Atji. Eq. |3.25| indicates that there is 
a preference for a positive total contribution at the 1-a level. In the CMSSM, At is 
typically negative due to RGE effects. Eq. shows that the sign of the non-SM 
contribution to the muon anomalous magnetic moment depends upon the sign of 
IxMi and /iM2. In the CMSSM, Mi and M2 are positive and so the combination of 
the [g — 2), and BR{B — > ^^7) constraints implies a preference for a definite sign 



of /i. Bayesian analyses [113], EM demonstrated that the current statistical evidence 



for yU > in the CMSSM is weaker than many may pre-suppose. 
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Figure 6: Posterior PDF of some of the ohservahles used to constrain the pMSSM. Solid 
(broken) lines represents plots for the linear (log) prior case. 



We have not used the freedom to re-define the phases of the fields and make M2 
positive for example, and so negative M2 appears in our fits. As such, in the pMSSM 
both At and Mi, M2 may take either sign and so the preference for yU > is broken 
and we may expect that the dominant contributions to the observables do not prefer 
either sign. On the other hand, there may be some residual dependence from the 
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sub-dominant contributions, as well as sub- dominant radiative corrections to other 
observables. Thus it is still important to check the relative probabilities for sign{^). 
For the pMSSM with a linear prior measure we find the following probability for the 
two signs of /i: 

P{fi > 0) = 0.40 ± 0.01 and P{fx < 0) = 0.60 ± 0.01. (4.5) 

Fig. ^ shows the posterior PDF marginalised on the M2 — fi plane. The figure shows 
that, although opposite signs for /i are allowed, it is constrained to have the same 
sign as M2, especially for log priors where {g — 2)^ is well fit. For linear priors, 
the large volume of parameter space leading to heavy sparticle masses mean that 
this tendency is reduced and there is a small amount of probability that /iM2 < 0, 
predicting a small negative 5a^. Fig. |^ clearly shows a symmetry of the fits when 
one simultaneously fiips the signs of M2 and /x, as should be the case due to phase 
re-definition freedom. 

I ^ ^ ^ ^ — , ^ 1 4000 F ' ' ' ' ' ' ' 




[GeV] [GeV] 

Figure 7: Marginalised 2D posterior PDF of M2 versus ^ for a (left) linear prior and 
(right) log prior. The dark contours show the 68% and 95% Bayesian credibility regions. 

The SM prediction of {g — 2)^ remains somewhat controversial. The hadronic 
contributions that are extracted from r and e~^e~ data disagree, and one obtains 
quite different 6a^ constraints depending upon which data set is used. To quantify 
the extent to which our mild preference for /i < depends on the {g — 2)^ observable, 
we made a separate sampling with all except the {g — 2)^ constraints. We found that 
/i < is still more probable with: 

P{fi > 0) = 0.46 ± 0.02 and P{fx < 0) = 0.54 ± 0.02. (4.6) 

Thus, {g — 2)^ contributes around 0.06 to the probability of yU < 0, the other observ- 
ables including BR[B — Xg'-y) contributing around 0.04. However, computing the 
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Bayesian evidence ratios in the two scenarios indicate that there are no conclusive 
evidence, based on Jeffrey's scale (see Tab. |l]), for one particular sign(/i) over the 
other. The odds and logarithm of the evidence ratios are summarised in Tab. 



Data Considered 


|loge 


Odds, Z+/Z_ 


Remark 


All constraints 

All minus [g — 2)^ constraints 


-0.41 ±0.04 
-0.18 ±0.04 


0.67 ±0.03 
0.84 ±0.04 


Inconclusive 
Inconclusive 



Table 8: The Bayesian evidence ratios for the pMSSM with /i > and ^ < with the linear 
priors measure. "All observables" refers to all the constraints discussed in Section [7^ . Z-^- 
and Z- represent the evidence for the hypothesis for the linear prior pMSSM with fi > 
and with n < respectively. A prior probability of 0.5 is assigned to each hypothesis. 



4.5 The gaugino code 

The LHC being a proton-proton scattering machine, is going to be producing a 
large number of strongly interacting particles. If TeV-scale SUSY is nature's choice 
of physics then the LHC machine would eventually be a gluino factory. The gluinos 
would cascade-decay down to the neutralino LSP. Thus the visible energy of the 
decay products of the gluino is determined by the gluino-LSP mass splitting. If the 
ratio of these two masses is large, there will be high visible energy and it should be 
easier to pick SUSY out from underneath backgrounds. The nature of the gluino and 
neutralino sparticles could also be important in discriminating the different models 



of SUSY. As recently emphasised [p.lO|| , patterns of gaugino masses can be used as 
an early discriminant of models of SUSY breaking at the LHC. The argument for 
this goes as follows. In, for instance, the CMSSM or minimal gauge-mediated SUSY- 
breaking, the boundary conditions of the values of the gauge coupling constants take 
the form gf : g2 : ^ 1 : 2 : 6 around the TeV scale. Here ga, a = 1,2, or 
3 respectively represent the electromagnetic, weak or strong interactions couplings. 
Unification of the gauge couplings at a GUT scale, gives a prediction on the pattern 
of gaugino mass terms since the ratio M^/g"^ does not run at one loop [ p.ll|| . If the 
neutralino is gaugino dominated, this translates into some gluino-to-neutralino mass 
ratio pattern. 

Already there are distinct ratios coming from viable SUSY breaking scenarios. 
For instance, the CMSSM (and AMSB) if it has a predominantly bino (and wino) 
LSP predicts rrig/m^o ^ 6 (and 9, respectively) [p.lO|| . Mirage mediation ||112|| with 
predominantly bino LSP and the large volume scenario (LVS ||113| , [114] , |115| , |116|| ) 



have a characteristic ratio less than 6 and between 3 to 4, respectively (with the 
LVS having the most compact gaugino mass pattern). Higgsino components of the 
neutralino LSP spoil a strict prediction of the mass ratio coming from gaugino mass 
ratios, however. By construction, the pMSSM set-up is the most generic and natural 
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Figure 8: pMSSM gluino-to-neutralino mass ratio for linear and log priors. 



approach for MSSM phenomenology. We show the pMSSM posterior PDF for the 
gluino-to-neutralino mass ratio in Fig. ^ which provides a clear discrimination be- 
tween the pMSSM and the other models. The figure ^ shows a sharp dependence of 
the gluino-to-neutralino mass ratio posterior PDF on the choice of prior. The linear 
prior predicts a compact gaugino mass ratio with*^ m^/m^o ^ 2.5 ± 1.0. For the log 
prior case a much broader distribution of mass ratios centred around m^/m^o ^ 10.0 
results from the fits. 

It is worthwhile to point out the source of the difference between the peak positions 
in Fig. ^ The gluino mass, m^, PDF is roughly the same for both priors. This implies 
that the numerator is fixed with respect to the two different priors. So the shift is 
solely due to the different properties of m^o which is much lighter for the log prior 
choice. 

4.6 Fine-tuning 

The main motivation of weak-scale SUSY is to solve the technical hierarchy prob- 
lem, explaining why the Higgs boson remains at the weak scale despite quantum 

^See for similar results from a pMSSM analysis done with a different set of observables than 
in this work. 
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corrections which are as large as the largest mass scale in the theory (e.g. the Planck 
scale). In order for softly broken SUSY to still provide a resolution of the technical 
hierarchy problem, the SUSY breaking terms should not be much larger than the 
TeV scale, otherwise an a priori un-natural cancellations between radiative correc- 
tions are required in order to keep the Higgs boson mass low. Direct SUSY search 
limits imply lower bounds on sparticle masses, which already start to imply that the 
MSSM parameters must cancel somewhat in order to separate the electroweak and 
SUSY breaking scales. This is termed the 'little hierarchy problem' [|1 1 7| , |118|| , (see a 
recent discussion in ||119|| ). We wish to quantify the necessary amount of fine-tuning 
in the pMSSM parameters needed to make the set-up consistent with the imposed 
sparticle mass bounds. 



We follow the approach in |p.20|| , quantifying the amount of fine-tuning in the 
Z-boson mass prediction coming from Higgs potential minimisation conditions. We 
consider this as a measure of fine-tuning in the pMSSM. The tree-level Z-boson mass 
is given by 
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The amount of fine-tuning is quantified by considering the sensitivity of niz to a 
variation of a parameter ^ 
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where ^ = m?^^,m\ and ^ are the relevant parameters in the pMSSM. Assuming 



tan /3 > 1, from Eqs. [4.7| [4.8| and [4.9| , one derives: 
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A(/i), A(m3), A(m|^J and ^{rri^jj^) are added in quadrature to obtain an over-all 
measure of fine-tuning, A^: 



Ar = ^A(^)2 + A(mi)2 + A(m2,J2 + A(m|,J2. 
Values of A^^ far greater than unity indicate large fine-tuning. 



(4.11) 
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Figure 9: Fine-tuning PDFs in the pMSSM. 



The posterior PDF for the amount of fine-tuning in the pMSSM is shown in Fig. |[ 
The logarithmic prior scenario have lower At than in the linear prior. This is not 
surprising since the SUSY breaking terms are much reduced in the former scenarios 
than in the latter. We see from the figure that the fine-tuning is most likely low at 
around At ~ 20 — 30, but there is a tail extending beyond At = 100. One could 



use a prior of 1/At, to encode a belief in less fine-tuned points in our global fits [0. 
Alternatively, one could place a cut on At, but the value of such a cut is of course 
subjective. 82% of the high posterior PDF points, around 4.0 x 10^ of samples, have 
At > 10, so a hard cut placed at 10 would have a drastic effect on the fits. Here we 
decline to change the prior or place cuts, since we are content with observing that, 
for most of the probability mass in the fits, it is not unacceptably large. For the 
highest likelihood models the fine-tuning is reduced from its average. For example, 
for the good-fit point in the linear prior sample. At = 24, whereas for the good-fit 
point in the log prior sample. At = 27. 

Notice that in general the amount of fine-tuning we find is small compared to pre- 
vious studies of the MSSM that start at a high scale, running down to the TeV scale 



using the RG equations. The reasons for this include, as explained in Refs. [122, 118 



that the amount of fine-tuning is a function of the cut-off scale and tends to decrease 
with this scale because the interval of RG running of the soft parameters induces 
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EWSB at tree level and the cross-talk (through RG running) between parameters 
in the Higgs sector and those from the squarks, gluinos, etc sectors is drastically 
reduced such that the latter parameters can be much heavier than mz without dis- 
turbing the naturalness of the electroweak scale^. Ref. | 122|| found some approximate 



semi-numerical solutions of the RGEs for the case that the boundary conditions on 
SUSY breaking parameters are set at the GUT scale ~ 10^^ GeV. The dominant 
term in /S.^{Mgut) typically comes from cross-talk with the GUT scale gluino mass 
M^{Mgut): 

Ai^iiMcuT)) = a(tan/3) ^'^ff^^^\ (4.12) 

where the authors determined the coefficient a(tan/3) numerically: a(2.5) = 24 and 
a(10) = 12, for example. This is to be contrasted with the pMSSM in Eq. [4.10 



where the terms in A(/x) are set at the SUSY breaking scale and are ~ 0(1), and do 
not involve the gluino mass, upon which there are strong empirical lower bounds. 

5. Neutralino Dark Matter 

Assuming R-parity conservation the lightest supersymmetric particle (LSP) may 
be a good dark matter (DM) candidate. In Ref. [^, a pMSSM study of the ability 



of SUSY measurements at future colliders to constrain dark matter properties was 
considered. We assume that the neutralino LSP constitutes the DM in the universe. 
The DM relic density then depends upon the LSP mass and, through its composi- 
tion in terms of gauginos and Higgsino, its interactions. 2-D marginalised posterior 
PDFs showing preferred regions in the relic density versus LSP mass are shown in 



Fig. [ly. There is a mild positive correlation between the preferred mass and the 
dark matter relic density QcDAih"^ for linear priors which is not evident for the log 
priors. It is clear that the LSP mass is not well constrained by current data, since 
it is highly prior dependent. The nature of the neutralino LSP in the pMSSM is ad- 
dressed in Section |5.1| . The mass difference between the LSP and the next to lightest 
supersymmetric particle (NLSP) is important because if it is small, the LSP may 
efficiently co-annihilate, in the early universe, with the NLSP, significantly reducing 
the relic density. The different possible NLSPs in the pMSSM with the correspond- 



ing posterior probabilities are listed in Section |5.2| . The dominant (co-)annihilation 



channels are presented in Section |5.2| . We present the prospects of direct dark matter 
detection in Section Ol 



5.1 Neutralino dark matter composition 

The nature of the neutralino LSP determines the (dominant) processes by which 
it (co-) annihilates into SM particles and therefore affects its present number density. 



^We thank J. R. Espinosa for interesting conversations on this point and on fine-tuning in general. 
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Figure 10: Marginalised 2D posterior PDFs of the neutralino dark matter relic density 
against the neutralino mass. Linear priors are on the left hand side and log priors on the 
right. The dark contours show the 68% and 95% Bayesian credibility regions. Note that 
the feature shown towards the left hand side of the linear prior plot is due to the peak in 
the neutralino mass at around 200 GeV in Fig. |^ which here just makes it into the 68% 
region. 



This was illustrated using the randomly scanned pMSSM samples in |^ where it was 
shown that the nature of the neutralino LSP depends upon whether one assumes 
that the LSP makes up all or only some of the DM relic density. The neutralino 



mass matrix is given by ^ip'^'^ M^ip^ + H.c. where ip 
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and. 



(5.1) 



Cx = cosx and Sx = sin a;. The neutralino mass eigenstates are Xi = where N is 

a unitary transformation that diagonalises M^. The LSP neutralino mass eigenstate 
is therefore a mixture of bino, wino and Higgsino: 



X 



Nub + Ni2W^ + Ni^H^ + Ni^Hl 



(5.2) 



Different regions of parameter space give different neutralino LSP compositions. 
When Ml -C min{M2, |/i|), Nu ~ 1 and the LSP is dominantly bino. Bino LSPs give 
a relic density that is too high for most of the parameter space unless some specific 
mechanism (such as efficient co-annihilations or annihilations through a resonance) 
is working. When M2 < min{Mi, |yu|), A^i2 dominates, i.e. the LSP is dominantly 
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wino and is quasi mass degenerate with the hghtest chargino. This leads to strong co- 
annihilations between the LSP and the chargino, and typically the relic density much 
smaller than the WMAP constraint for wino LSPs. For < mm(Mi, M2), N13 and 
Ni4^ are of order one and the LSP is dominantly Higgsino and there may be efficient 
annihilations into top and weak gauge boson pairs. In the Higgsino-dominated LSP 
scenario, Xi; X2; xf are almost all mass degenerate and are Higgsino- like. Of course, 
there exist mixed cases which include several of these limiting behaviours. 




iog,„[i-z^] 



Figure 11: pMSSM neutralino gaugino-Higgsino admixture fractions. Higgsino domina- 
tion is at the right-hand side of the plot and gaugino domination is to the left-hand side. 



The gaugino/Higgsino mixture PDF of the LSP is shown in Fig. |TT| constructed 
from the fraction 

Z,= \Nu\'+\N22\' 



following [ [123|| . 1 — Zg is unity if the neutralino LSP is fully Higgsino-like and zero if 



fully gaugino-like. The plot shows that LSP is mostly Higgsino in the linear prior case. 



similar to the non- universal higgs mass scenario ||124|| , and mostly gaugino for the log 
prior scenario. Thus, current data do not unambiguously constrain the LSP content. 
Referring to Tab. ^ we see that the good-fit point from the linear prior sample 
has a mixed wino-Higgsino LSP (more precisely, the point has {Ni^l ~ |-^i4| = 0.7 
and = 0.15). The log prior sample good-fit point has a bino dominated LSP 
(1^11 1 = 0.998), but there are several light sparticles, allowing sufficient annihilation. 

5.2 (Co-) Annihilations 

At early times of the universe the LSP is in thermal equilibrium with other particles 
and, ignoring for now co- annihilations, its number density evolution is governed by 
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the Boltzmann equation 



_^ = -3i7n^o- <av> {nl. - (ng)^}. (5.3) 

Here H is the Hubble expansion rate of the universe, n^o is the number density and 
< av > is the thermally averaged annihilation cross sections of the neutralino LSP. 
V is the relative velocity of the annihilating pair. The LSP annihilation rate is given 
by r^o =< av > n^o. At a freeze-out temperature Tj, the neutralino decouples, 
r^o = H{Tf). Substituting H{Tj) into Eq. |5]^ predicts that the LSP rehc density is 
inversely proportional to the thermally averaged annihilation cross section, < av >. 
This means that for the LSP relic abundance today to be in the WMAP-5 range 
Eq. |3.38| there must be a significant number of annihilations of the neutralino LSP at 



earlier times. The possible processes are mostly two-particle final states which could 
be fermion anti-fermion pair, combinations of the weak gauge-bosons or 
combinations of the Higgs bosons {hP , , , H^) (see e.g. ||125|| ). The discussion 
becomes much more involved once co-annihilation processes are taken into account, 
since coupled Boltzmann equations are required for each relevant SUSY particle 
species. 

Co-annihilation processes dominate in parameter space points where next-to-lightest 
supersymmetric particle (NLSP) that are almost mass degenerate with the LSP. At 
such points the neutralino abundance also depends strongly on the annihilations of 
the NLSPs ||126|, p.27|| and the number densities of the LSPs and NLSPs are coupled. 



A review of different co-annihilation studies was presented in Here we present 
and analyse the outcome of the pMSSM annihilation and co-annihilation results for 
our two different prior measures. We shall only discuss processes that contribute 1% 
or more of the annihilation cross section. 

In Tab. ^we give a list of possible NLSPs and corresponding posterior probabilities 
for each. The probabilities indicates that neutralino-chargino co-annihilations are 
most likely to be dominant in the pMSSM with a linear prior measure. For the 
log prior measure, neutralino-slepton co-annihilations are the most probable. We 
illustrate the abundance and complexity of the annihilation channels by presenting 



the dominant ones at the good fit points in Tab. 10. The dominant channels for the 
linear prior sample good-fit point are direct neutralino-chargino co-annihilation and 
neutralino annihilation via chargino exchange into Z- and W-hoson pairs. For the log 
prior measure good- fit point, the dominating channels are neutralino co- annihilations 
with various sleptons. Many different processes contribute at the percent level. We 
present the identities of the dominant channel, along with its posterior probability in 



Tab. IL The most likely channel is neutralino-chargino co-annihilation for the linear 



prior and neutralino annihilation for the log prior. From the large prior dependence 
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in the results, we deduce that current data are not sufficient to constrain the dark 
matter annihilation properties of the LSP. 
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Table 9: pMSSM NLSP identity probabilities for linear and log priors. 



5.3 Direct detection 



128 



Many different experiments search for the nature of dark matter (e.g. see 
and references therein). Indirect detection experiments are designed to observe the 
annihilation products of dark matter particles. We do not address indirect detection 
here and save it for future consideration. Here, we consider direct detection experi- 



ments such as XENON pg, CDMS |[T30|, |m|], ZEPLIN [|132| [T3|, Edelweiss |34 



CRESST |T35[, WARP |T36|, g or COUPP [g. Such experiments are designed 
to observe the elastic scattering of dark matter particles with nuclei. The LSP may 
interact with quarks in target nuclei via t-channel CP-even Higgs exchange or s- 
channel squark exchange and with gluons via squark loop contributions. DM direct 
detection rates also depend on the local neighbourhood DM density and velocity 
distribution. The density, which is estimated to lie between 4 x 10~^^ g/cm~^ and 
13 X 10~'^^g/cm~^ (0.22 — 0.73 GeV/cm^), is inferred by fitting observations to mod- 
els of galactic halo ||139| , |140|| . The velocity is expected to be around 230 ± 20 km/s 



141 



The elastic scattering cross section is partitioned into spin-dependent and spin- 
independent components. The spin independent part is currently the most con- 
straining, and we concentrate on it. It is proportional to the square of the target 
nucleus atomic number, A^. This enhancement is because the dark matter wave- 
length is of same order as the size of a nucleus and hence the scattering amplitudes 
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Table 10: (Co-) annihilation channels at the good fit points. We display the percentage 
contribution to the annihilation cross section for each channel. Channels which contribute 
less than 1% to the annihilation cross sections are not listed. 
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Table 11: Posterior probabilities for dominant annihilation and co- annihilation channels. 



on individual nucleons add coherently. There is one experimental claim in direct 
detection experiments of a signal in the annual modulation rate ||142|| . This result 
has not been confirmed by other experiments and would be incompatible with a neu- 
tralino LSP candidate, so we do not use it to constrain the pMSSM. Aside from this, 
no positive signal has been detected to date in dark matter detection experiments. 
A positive signal would constrain SUSY parameter space if one assumed a particular 
local neighbourhood DM density and velocity distribution. 

The spin-independent neutralino-nucleus elastic scattering cross section is given 
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by 



a 



:[Zfp + {A-Z)U 



(5.4) 



where rriT is the mass of the target nucleus and Z and A are the atomic number and 
atomic mass of the nucleus, respectively, fp and /„ are the neutralino's couplings to 



protons and neutrons, given by ||143 

E ft'" 
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where are the neutralino-quark couplings [|Tg |n6|, |Ti7[ ^ and 

denote the quark content of the nucleon. They have been experimentally bounded 



to be: /. 



(p) 



0.020 ± 0.004, /, 



(p) 



0.026 ± 0.005, /. 



(p) 



0.118 ± 0.062, /, 



(n) 



0.014 ± 0.003, f^^ ^ 0.036 ± 0.008 and /j,"^ ^ 0.118 ± 0.062 m, |l50l The 



first term in Eq. |5.5| corresponds to interactions with the quarks in the target, which 
can occur through either t-channel CP-even Higgs exchange, or s-channel squark 
exchange. The second term corresponds to interactions with gluons in the target 
through a quark/squark loop. f^j^Q is given by 1 — f!^J — f^^ — f^^ ~ 0.84, and 



analogously, ~ 0.83. 




Figure 12: Posterior PDF of the neutralino-proton spin-independent scattering cross sec- 
tion for the pMSSM with linear (left) and log (right) prior measures. The CDMS 90% 
confidence level upper hound is also shown, assuming a local DM density of 0.3 GeV/cnr? . 
The dark contours show the 68% and 95% Bayesian credibility regions. 

The direct detection constraints from the cryogenic cold dark matter search (CDMS) 
experiments on the pMSSM is shown in Fig. The large prior dependence of the 
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results indicates that current data are insufficient to constrain the direct detection 
cross sections. One can say that there is clearly a wide allowed range for the direct 
detection cross sections. 



5.4 Relaxing the purely LSP dark matter assumption 



Gaussian 
upper bound 




0.2 



^CDM h 



Figure 13: Depiction of the likelihood constraint on the predicted value of Qdm^^ due to 
lightest neutralinos, compared to a simple Gaussian with WMAP5 central value and a la 
uncertainty of 0. 02 used in the rest of the paper. 



The analysis presented above was done assuming that the neutralino LSP is the 
only source of dark matter. It is known that the LSP relic density is sensitive to 
the assumed cosmology. For example, Big Bang Nucleosynthesis (BBN) expansion 
rates can enhance the calculated relic density without affecting other important cos- 
mological quantities ||152[ . The inclusion of right-handed neutrinos could change the 
relic density prediction, see for instance [|153|] . One could also allow for additional 
non-neutralino dark matter components. In order to see the potential effect of such 
model changes, we relax the constraint from the DM relic density to the case where 
only QcDAih'^ predictions larger than the central values are penalised according to 
the likelihoods: 



where /x = 0.1143 is the experimental central value quoted above and s is an inflated 
error on relic density that includes theoretical uncertainties in its prediction. 



We have made an independent run using Eq. p.6| , i.e. relaxing the purely LSP DM 
assumption, i.e. implicitly assuming some other component of dark matter. We wish 
to examine the amount of DM that comes from the LSP. These runs were performed 
in Ref. ||3^, where the relevant constraints can be found. Linear priors were used 
on the parameters, which had a 2 TeV upper bound. We find in Fig. |l^ that the 
preferred relic density is low compared with the purely LSP DM assumption: around 
QcDMh"^ = 10^^ — 10~^. Thus once one allows for an additional component of DM to 
the LSP, the model prefers the additional component to dominate the relic density 

10 



with full WMAP5 Gaussian constraint 
with only WMAP5 upper bound constraint 



log a.h^ 

10 X, 



Figure 14: Neutralino relic density assuming WMAP5 as a Gaussian likelihood constraint 
or as an upper bound. For this plot only, a 2 TeV range linear pMSSM with other param- 
eters as in Ref. was taken. 



6. Conclusions and Outlook 

We have presented the first statistically convergent global fit of the pMSSM model 
with its 25 independent continuous parameters plus a discrete parameter, sign^jj) = 
±, to the dark matter relic density, indirect observables and direct sparticle search 
constraints. We have used the entire set of relevant electroweak precision observables 
and 5-physics data, as well as the anomalous magnetic moment of the muon as 
indirect observables. The evidence for the linear and log prior measure pMSSM in 
light of the data is log^ Z = 63.211 ± 0.033 and log, Z = 65.043 ± 0.042 respectively. 



-^°Some of us hope to return to this issue in future work. We thank Bryan Webber for suggesting 
this comparison. 
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We have presented a good-fit point in the pMSSM parameter-space which is not 
unique since the minimum is degenerate. All the different parameter-space points 
with maximal likelihood are equally interesting for phenomenological studies. 

This work constitutes a technical demonstration that statistically convergent global 
fits in high dimensions involving curving degeneracies and several modes are now fea- 
sible. This feasibility is due to new sampling algorithms and improvements in the 
speed of computation and access to it. It allows more complete phenomenological 
studies of multi-parameter models beyond the standard model that could not have 
been completed in the past. In particular, the set-up and techniques employed here 
provide an unbiased approach to MSSM phenomenology - independent of the under- 
lying theory, the mechanism to break SUSY or its mediation - hence it could lead 
to more robust SUSY phenomenological studies and guides for LHC SUSY searches 
and for dark matter search experiments. As expected the results of this exercise 
differs significantly from those of the more studied CMSSM/mSUGRA and thus the 
pMSSM parameter space provides a much richer and appropriate arena for LHC 
studies of the MSSM. 

For the analysis we consider prior measures flat in either the parameters (a lin- 
ear prior) or flat in the logarithm of the parameters (a log prior) in order to check 
robustness of inferences. Given the large number of (the pMSSM) parameters com- 
pared with the (weak constraining power of the) data available at the moment it is 
very interesting that indeed there are some prior-independent results or inferences. 
The lightest CP-even Higgs boson mass and the stop masses fall into this category of 
robust and prior independent results. The other sparticle masses and all dark matter 
properties exhibit significant prior dependence and require more direct and precise 
data (or more constrained models) to make their prediction robust. We emphasise 
that prior dependence is present and it is a positive feature of Bayesian methods, 
since its absence signals when there is enough data to make the fits robust. This 
can be used as a guidehne for future studies of the experimental implications of 
the MSSM. In particular, if SUSY is discovered via sparticle production at collider 
and many sparticle properties are precisely measured, it will be possible to use the 
techniques show-cased in the present work to extract pMSSM information. Such 
information can then be checked for consistency with more constrained models. 

We now contrast our methods with the recent random points scan analysis of the 
pMSSM in Ref. and a similar earlier work in Ref. [Q. These works perform a 
pMSSM random parameter scan to find points which pass direct search and dark 
matter constraints while being within 2 a of the central values of some indirect 
constraints. All such points are considered on an equal footing, and as such (as 
emphasised by the authors) are not a statistical global fit. We, on the other hand, 
allow a trade-off between different observables in a statistically correct fashion; one 
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may tolerate a moderately bad fit in one observable if it fits the other observables 
particularly well. Our use of Gaussian distributions for the likelihoods (instead of 
the 2(7 top-hat function in M, R3|) is justified by the central limit theorem and the 



maximum entropy principle ||154|| . Interesting points with LHC phenomenology not 
covered by previous studies of constrained models were found in Ref. [Q, which was 
the main aim of the approach (a few thousand points passed the constraints, out of 
~ 10^ scanned). We aim to perform a complete and statistically convergent global fit 
of the pMSSM. To achieve this we take advantage of the power of the MultiNest 
algorithm, which provides samples in moderately high dimensional parameter spaces 
(with curving degeneracies and different modes) much more efficiently than in ran- 
dom/grid parameter scans. In ||3^ there was more emphasis on direct search limits. 



which are more sophisticated than the ones employed in the present paper. Since 
the sparticle masses implied by our fits are large, our results are insensitive to the 
exact form of the direct search limits. The density of points in Ref. |]3^ also shows 
prior dependence although the results were not interpreted statistically (and thus 
they were not Bayesian nor frequentist). Another major difference to our approach 
is that in Ref. [Q, the WMAP constraint is used only as an upper bound (making 
viable points much easier to find) and so the existence of another non-MSSM dark 
matter particle is assumed. This changes the character of the points: allowing MSSM 
points which predict approximately zero LSP relic density means that a large num- 
ber of sampled points have a wino dominated LSP. We, however, assume in most of 



our analysis that the neutralino constitutes all of the dark matter. In Section ( |5.4| ) 
we presented an independent run made using the WMAP constraint only as a lower 
bound and our results (as expected) agree with those of [^, ^] in the sense that, in 
that case, the LSP contribution tends to be only a small fraction of the total dark 
matter. 

There are many directions in which this research could be extended. For each of 
the tens of thousands of preferred points in our sample, a detailed calculation of LHC 
observables such as inclusive counts of opposite sign dilepton and trilepton events, 
could be made using standard event generators and detector simulators, as has been 
done for more constrained models such as the CMSSM ||155||. This would provide a 



portrait of the signature space that may eventually be useful in direct SUSY searches. 
On a simpler level, one could compute relative probabilities of various sparticle mass 
hierarchies. The indirect dark matter detection prospects could also be evaluated, 
although with current data, they are likely to be highly prior dependent. The impact 
of the inclusion of the fine-tuning into the prior could be analysed. Assuming a 
particular parameter point, the impact of LHC SUSY measurements on our fits 
could be evaluated and an estimate of how much luminosity would be required in 
order to make inferences approximately prior independent^^. In this case, model 
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A Bayesian approach has recently been used to amehorate the LHC inverse problem in the 
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comparison between more constrained models and the pMSSM could be informative. 
One could determine, for a given LHC luminosity, which models could be made prior 
independent. The extension of the analysis to the full 124 MSSM parameter space 
may be still out of reach at the moment. Algorithms improving the MultiNest 
algorithm may be required before attempting it. Adding a small number of minimal 
flavor violating parameters could however be feasible. An extension of this work to 
include reasonable generalisations of the minimal flavour violation scenario adding a 
few extra parameters should be possible. Also, including R-parity violation to the 
pMSSM as well as a phenomenological NMSSM are within reach of the techniques 
we used here. Some of us hope to return to these issues in future work. 



A. Nested sampling and the MultiNest algorithm 

For scanning parameter spaces of large dimensionality we have to use efficient 
modern approaches for sampling the posterior. In such problems, interesting param- 
eter regions are often tiny in some directions and many directions are orthogonal 
to ones along which the likelihood is degenerately high. In this Section we present 
the procedure, in the context of the pMSSM, for the Monte Carlo technique called 
nested sampling developed by John Skilling and implemented in MultiNest. 
It is a general method for evaluating the integral Eq. |2.2| from which representative 
samples from the posterior distribution Eq. are obtained as by-product. The 
method differs from the traditional approach to inference dating back to Metropolis 
et. al. (1953) where the emphasis is more on evaluating the posterior density than 
in calculating the evidence. Skilling's method goes as follows. Exploring the 25- 
dimensional co-ordinate to evaluate the evidence integral is impractical. Instead, 
the prior mass dX = 7r(0) d@ can be used directly to convert the 25-dimensional into 
a 1-dimensional integral over a unit interval. Let X{L) be the prior mass enclosed 
within the likelihood contour, L{@) = L in the parameter space. That is, 

X{L) = [ 7r(e)rf25e. (A.l) 

Jl{&)>l 

As L increases from zero to infinity the enclosed prior mass decreases from X{0) = 1 
to X(oo) = 0. The inverse function L{X) = L is the contour value (a likelihood 
value) such that the volume enclosed is X (see Fig. |l^ for an illustration). Eq. [A.l| 



and the definition of its inverse implies that the evidence Eq. ^]2| can be expressed 
as ^ 

Z = j L{@)7i{@) d^^e = L{X)dX. (A.2) 
Given the likelihood values Li = L{Xi) at a sequence of m points < Xm < . . . < 
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Figure 15: Likelihood contours over a toy two-parameter space showing the enclosed vol- 
ume mapped to corresponding prior mass. A prior mass x = ^ is mapped to the likelihood 
contour that encloses i of the prior volume. Note the nested nature of the contour lines. 



Figures are from [ 31 /. 



X2 < Xi < Xo = 1 the evidence is estimated as a weighted sum, 

in 

Z = Y,L^w,, (A.3) 
1=1 

where for the trapezoidal rule Wi = ^(Xj_i — Xj+i). 

A.l Evidence evaluation 

The nested sampling procedure for evaluating the evidence starts with the ac- 
cumulation of N points uniformly drawn from the prior, the initialisation of the 
evidence, Z = 0, and the initialisation of the prior volume, Xq = 1. The number, 
N, of "live" points, 0i,...,0Ar is preserved throughout the procedure and every 
point is associated with its corresponding likelihood value: L(0i), . . . , L{@iy). Each 
step i = 1, 2, 3, . . . over the iterations is associated with the lowest likelihood Li (or 
the largest prior mass, Xj) that defines the contour line (or shell) L{@) = Li over 
parameter space. For moving from the {i — l)th to the iih iteration a new point is 
drawn from the set of points uniformly distributed in (0,Xj_i), the parameter space 
region with likelihoods L > L(Xj) = Lj. This is illustrated in Fig. 0. The new point 
replaces the one with lowest likelihood. Xj is set to Xj = exp{—i/N), the weight 
Wi to |(Xj_i — Xj+i) and the evidence Z incremented by LiWi. This procedure is 
repeated for the subsequent iterations. 
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Figure 16: Figure on the left shows N = 8 live points uniformly distributed in parameter 
space (or prior volume space (0, 1)) and sorted according to corresponding likelihood values. 
On the right is a picture illustrating the sampling of a new point (big purple) dot from the 
live points uniformly distributed in (0,a;i). Figures are from 



The prior volume shrinkage ratios ti = Xj/Xj_i are distributed according to 
Pr{ti) = Nt^~^ in (0, 1) where ti is the largest of random numbers uniformly 
distributed in (0,1). Sampling over t represents a geometrical exploration of the 
parameter space. The mean and standard deviation of t are 

E{\ogt) = and a[\ogt] = ^ (A.4) 

respectively. This justifies the assignment Xi = exp{—i/N) since each draw of logt 
is independent and after i iterations of the sampling procedure the prior volume will 
shrink down according to 

hgXi^ -{i±Vi)/N. (A.5) 

A. 2 Stopping criterion 

The nested sampling procedure is terminated after a preset number of the iterations 
(as described in Subsection |A.1|) or when the largest likelihood taken over the whole 
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currently (at the instance of check for termination, say the jth iteration) available 
prior mass would not increase the evidence value by more than some preset fraction 
/ (we use 0.5 in log-evidence). That is, the procedure is terminated if 

max(L(0i), . . . , L{@n))Xj < fZj. 

The integration, Z, is dominated around the region X ^ , wherever the bulk 
of the posterior mass is to be found. Here 

H = information = - j \og{dX/dP) dP^Yl log (^^^ (A.6) 

and dX/ dP is the compression ratio representing the fraction of the prior mass that 
contains the bulk of the posterior. dP = p{@) d@ = Z^^L{@)7r{@),d@. Recalling 
that Xi ^ e~*/^ we expect the integration procedure to take NH ± a/ NH steps 
(iterations) before reaching covering the bulk of the posterior. Hence another ter- 
mination condition could be to continue iterating until the count i is significantly 
greater than NH. 



The uncertainty in X translates to a geometrical uncertainty factor exp{±>y H /N) 
in the weights Wi of the dominating iterates. This in turn gives the uncertainty in 
the evidence via Eq. ^ as devilogZ) ^ ^H/N so that 



\ogZ = \og[J2L^wA±J^. (A.7) 



■t=i 



A. 3 Posterior inferences 

The posterior distribution p{@) is simply the prior distribution weighed by the 
likelihood. This can be trivially extracted from the evidence calculation since the 
set of sampled points 0i, . . . , ©at is already a posterior representative provided it is 
assigned the appropriate importance weight and normalised by the evidence, Z, to 
produce probability density with unit total. That is at the ith iteration the posterior 
probability density is 

LiWi . , 

Pi = (A.8) 

These are generated from the sequences of discarded points (the points with the 
lowest likelihood value at each iteration) during the sampling procedure. From these 
posterior sequence properties such as the mean /i and standard deviation a of some 
Q{@) are easily computable: 

/i = $^P^Q(e,) and fx' + a' = J^PiQi^if- (^-9) 

i i 

Equally-weighed samples selected proportionally to Pi can be used to construct 
marginalised posterior distributions in 0. 
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(b) 




Figure 17: Illustrations of the ellipsoidal decompositions performed by MultiNest. The 
points given as input are overlaid on the resulting ellipsoids. 1000 points were sampled 
uniformly from: (a) two non-intersecting ellipsoids; and (b) a torus. 

For completeness, it is worth mentioning that there are alternative methods for 
evaluating the evidence with other advanced MCMC algorithms like thermodynamic 
integration and it is not clear yet which method is best for high dimensional prob- 
lems: dimensions greater than 10. However, for this paper we implement the 
nested sampling algorithm for our purpose using the MultiNest code which 
has additional quality of being efficient in sampling multi-modal posteriors exhibiting 
curving degeneracies (see a summarised account in Subsection |A.4| ). 



A. 4 MultiNest 

The most challenging task in implementing the nested sampling algorithm is draw- 
ing samples from the prior within the hard constraint L > Lj at each iteration i. Em- 
ploying a naive approach that draws blindly from the prior would result in a steady 
decrease in the acceptance rate of new samples with decreasing prior volume (and 
increasing likelihood). MultiNest algorithm ^ tackles this problem through 



an ellipsoidal rejection sampling scheme by enclosing the live point set into a set of 
(possibly overlapping) ellipsoids and a new point is then drawn uniformly from the 
region enclosed by these ellipsoids. The number of points in an individual ellipsoid 
and the total number of ellipsoids is decided by a an 'expectation-maximisation' 
algorithm so that the total sampling volume, which is equal to the sum of volumes 
of the ellipsoids, is minimised. This allows maximum flexibility and efficiency by 



"'^^We thank David Mackay for interesting discussions about this. See for example 
http : / / www . inference . phy. cam. ac . uk / mackay /present ations /nestedOG/ 
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breaking up a mode resembling a Gaussian into relatively fewer number of ellip- 
soids, and if the posterior mode possesses a pronounced curving degeneracy so that 
it more closely resembles a (multi-dimensional) 'banana' then it is broken into a 
relatively large number of small 'overlapping' ellipsoids (see Fig. |1^. With enough 
live points, this approach allows the detection of all the modes simultaneously re- 
sulting in typically two orders of magnitude improvement in efficiency and accuracy 
over standard methods for inference problems in cosmology and particle physics phe- 
nomenology (see e.g. ||157| , |158| , |34| , ^ |159|] ). The MultiNest procedure as applied 



to our pMSSM fits is summarised by the fiow charts in Fig. Ill 
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Initialisation: N=4000 live points 
and corresponding likelihoods 



Z = 0;X[0]=1 
counter = 



counter = counter + 1 




pick lowest likelihood object L[i] 
get corresponding weight w[i] 
update Z = Z + L[i] w[i] 

posterior p[i]=L[i] w[i]/Z 

[ 

update and save file with 
p[i], L[i]. and the parameter points 



replace object L[i] with new one 
Flowchart 3 



Main flowchart 




compute likelihood 
keep point and likelihood 




Draw a 26-D point from L( 9 ) > L 



Pass parameter to HEPsoftwares 




Replace point 




Figure 18: Flow charts summarising the sampling procedure: Z refers to the Bayesian 
evidence, Eqs. (2^, (AA^), (A_^); Xi is the prior mass, Eq. A.l ; and pi is the posterior 
probability, Eq. A.^ . HEP software refers to the different computer programs described in 
Section ^3). 
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